# Get solutions to the gauged GL vortices in (1+1) dimensions?

Posted 10 months ago
1145 Views
|
8 Replies
|
5 Total Likes
|
 I am trying to numerically solve the following differential equations for the profiles of the Gauged GL Vortices: eqn1 = y''[x] + y'[x]/x - ((k - z[x])^2* y[x])/x^2 + ((1 - y[x]^2)*y[x])/2 == 0; eqn2 = z''[x] - z'[x]/x + (k - z[x])*y[x]^2 == 0; with the boundary conditions: inf = 1000; bc = {y[0] == 0, y[inf] == 1, z[0] == 0, z[inf] == 1}; I have tried a finite difference method but this seems to give large errors even for k = 1 so I am wondering if there is a better numerical method to yield a more accurate solution?Below is my current function which I employ to solve these coupled differential equations: Clear[fdd, pdetoode, tooderule, sollst] fdd[{}, grid_, value_, order_] := value; fdd[a__] := NDSolveFiniteDifferenceDerivative@a; pdetoode[funcvalue_List, rest__] := pdetoode[(Alternatives @@ Head /@ funcvalue) @@ funcvalue[[1]], rest]; pdetoode[{func__}[var__], rest__] := pdetoode[Alternatives[func][var], rest]; pdetoode[rest__, grid_?VectorQ, o_Integer] := pdetoode[rest, {grid}, o]; pdetoode[func_[var__], time_, {grid : {__} ..}, o_Integer] := With[{pos = Position[{var}, time][[1, 1]]}, With[{bound = #[[{1, -1}]] & /@ {grid}, pat = Repeated[_, {pos - 1}], spacevar = Alternatives @@ Delete[{var}, pos]}, With[{coordtoindex = Function[coord, MapThread[ Piecewise[{{1, # === #2[[1]]}, {-1, # === #2[[-1]]}}, All] &, {coord, bound}]]}, tooderule@ Flatten@{((u : func) | Derivative[dx1 : pat, dt_, dx2___][(u : func)])[x1 : pat, t_, x2___] :> (Sow@coordtoindex@{x1, x2}; fdd[{dx1, dx2}, {grid}, Outer[Derivative[dt][u@##]@t &, grid], "DifferenceOrder" -> o]), inde : spacevar :> With[{i = Position[spacevar, inde][[1, 1]]}, Outer[Slot@i &, grid]]}]]]; tooderule[rule_][pde_List] := tooderule[rule] /@ pde; tooderule[rule_]@Equal[a_, b_] := Equal[tooderule[rule][a - b], 0] //. eqn : HoldPattern@Equal[_, _] :> Thread@eqn; tooderule[rule_][expr_] := #[[Sequence @@ #2[[1, 1]]]] & @@ Reap[expr /. rule] Clear@pdetoae; pdetoae[funcvalue_List, rest__] := pdetoae[(Alternatives @@ Head /@ funcvalue) @@ funcvalue[[1]], rest]; pdetoae[{func__}[var__], rest__] := pdetoae[Alternatives[func][var], rest]; pdetoae[func_[var__], rest__] := Module[{t}, Function[pde, #[ pde /. {Derivative[d__][u : func][inde__] :> Derivative[d, 0][u][inde, t], (u : func)[inde__] :> u[inde, t]}] /. (u : func)[i__][t] :> u[i]] &@ pdetoode[func[var, t], t, rest]]
8 Replies
Sort By:
Posted 10 months ago
 The following could be used as a starting point for solving the system numerically: Clear["Global*"]; k = 2.; inf = 1000.; xg = N@Range[0, inf, inf/12000]; xg1 = xg + 10^-50.; {eye, dX, d2X} = NDSolveFiniteDifferenceDerivative[#, xg, "DifferenceOrder" -> 6]["DifferentiationMatrix"] & /@ {0, 1, 2}; mu = -3*10^2.; mat = d2X + mu*eye; mat[[{1, -1}]] = eye[[{1, -1}]]; lfun = LinearSolve[mat]; ysol = zsol = Erf[xg/3]; Do[ rhs1 = (((k - zsol)^2)*ysol/xg1^2) - ((1 - ysol^2)*ysol)/2 - (dX.ysol)/xg1 + mu*ysol; rhs1[[{1, -1}]] = {0., 1.}; rhs2 = -((k - zsol)*ysol*ysol) + (dX.zsol)/xg1 + mu*zsol; rhs2[[{1, -1}]] = {0., 1.}; ysol1 = lfun[rhs1]; zsol1 = lfun[rhs2]; If[Max[Abs[ysol1]] >= 10 || Max[Abs[zsol1]] >= 10, Print["unstable.."]; Break[]]; err = Max[Norm[ysol1 - ysol], Norm[zsol1 - zsol]]; If[err <= 10^-7, Break[]]; If[Mod[iter, 1000] == 0, Print["err ->", err]]; ysol = ysol1; zsol = zsol1; , {iter, 1000000}] The parameter mu will have to be changed as the parameter k is changed. A larger value of k would require a larger value of mu. Also, for larger k-values, the convergence associated with these iterations is slow.
Posted 9 months ago
 Dear Paritosh Mokhasi, Is there any description/explanation of this method of solving the differential equations you could provide? For example a paper or link explaining it.
Posted 9 months ago
 The method I am using here is called Successive Over Relaxation (SOR). There are plenty of references that you can find on the web on this. My implementation is a very simple version. There are much more sophisticated methods in SOR that could be used.
Posted 10 months ago
 Unfortunately, the numerical solutions associated with your method are actually worse than my original finite difference method which I have shown above, rather than better.
Posted 10 months ago
Posted 10 months ago
 Clear[fdd, pdetoode, tooderule, sollst] fdd[{}, grid_, value_, order_] := value; fdd[a__] := NDSolveFiniteDifferenceDerivative@a; pdetoode[funcvalue_List, rest__] := pdetoode[(Alternatives @@ Head /@ funcvalue) @@ funcvalue[[1]], rest]; pdetoode[{func__}[var__], rest__] := pdetoode[Alternatives[func][var], rest]; pdetoode[rest__, grid_?VectorQ, o_Integer] := pdetoode[rest, {grid}, o]; pdetoode[func_[var__], time_, {grid : {__} ..}, o_Integer] := With[{pos = Position[{var}, time][[1, 1]]}, With[{bound = #[[{1, -1}]] & /@ {grid}, pat = Repeated[_, {pos - 1}], spacevar = Alternatives @@ Delete[{var}, pos]}, With[{coordtoindex = Function[coord, MapThread[ Piecewise[{{1, # === #2[[1]]}, {-1, # === #2[[-1]]}}, All] &, {coord, bound}]]}, tooderule@ Flatten@{((u : func) | Derivative[dx1 : pat, dt_, dx2___][(u : func)])[x1 : pat, t_, x2___] :> (Sow@coordtoindex@{x1, x2}; fdd[{dx1, dx2}, {grid}, Outer[Derivative[dt][u@##]@t &, grid], "DifferenceOrder" -> o]), inde : spacevar :> With[{i = Position[spacevar, inde][[1, 1]]}, Outer[Slot@i &, grid]]}]]]; tooderule[rule_][pde_List] := tooderule[rule] /@ pde; tooderule[rule_]@Equal[a_, b_] := Equal[tooderule[rule][a - b], 0] //. eqn : HoldPattern@Equal[_, _] :> Thread@eqn; tooderule[rule_][expr_] := #[[Sequence @@ #2[[1, 1]]]] & @@ Reap[expr /. rule] Clear@pdetoae; pdetoae[funcvalue_List, rest__] := pdetoae[(Alternatives @@ Head /@ funcvalue) @@ funcvalue[[1]], rest]; pdetoae[{func__}[var__], rest__] := pdetoae[Alternatives[func][var], rest]; pdetoae[func_[var__], rest__] := Module[{t}, Function[pde, #[ pde /. {Derivative[d__][u : func][inde__] :> Derivative[d, 0][u][inde, t], (u : func)[inde__] :> u[inde, t]}] /. (u : func)[i__][t] :> u[i]] &@ pdetoode[func[var, t], t, rest]] Clear@k; Clear@y; Clear@z; Clear@\[Phi]; inf = 1000; eqn1 = y''[x] + y'[x]/x - ((k - z[x])^2* y[x])/x^2 + ((1 - y[x]^2)* y[x])/2 == 0; eqn2 = z''[x] - z'[x]/x + (k - z[x])*y[x]^2 == 0; eqns = Join[{eqn1, eqn2}]; bc = {y[0] == 0, y[inf] == 1, z[0] == 0, z[inf] == k}; set = Join[eqns, bc]; points = 1000; domain = {0, inf}; grid = Array[# &, points, domain]; difforder = 4; (*Definition of pdetoode is above.*) ptoafunc = pdetoae[{y[x], z[x]}, grid, difforder]; delete = #[[2 ;; -2]] &; ae = delete /@ ptoafunc[Map[x^2 # &, eqns, {2}] // Simplify]; WorkingPrecision -> 30; sollst = Table[ ListInterpolation[#, grid] & /@ Partition[#, points] &@ With[{initial = 1}, FindRoot[{ae, bc}, Flatten[{{y@#, initial} & /@ grid, {z@#, initial} & /@ grid}, 1]]][[All, -1]], {k, 3}] Attachments: