1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
|
BeginPackage["Schofield`"]
$Assumptions = {θc > 0, θc > 1, gC[_] ∈ Reals, B > 0, γ > 0}
β[D_:2] := Piecewise[
{
{1/8, D == 2},
{0.326419, D == 3},
{1/2, D == 4},
{β, True}
}
]
δ[D_:2] := Piecewise[
{
{15, D == 2},
{4.78984, D == 3},
{3, D == 4},
{δ, True}
}
]
α[D_:2] := Piecewise[
{
{0, D == 2},
{0.11008, D == 3},
{0, D == 4},
{α, True}
}
]
Δ[D_:2] := β[D] δ[D]
f[θ_] := θ^2 - 1
g[gC_:gC, θc_:θc][n_][θ_] := (1 - (θ/θc)^2) Sum[gC[i] θ^(2i+1), {i, 0, n}]
I\[ScriptCapitalM]f[y_] := (1 + 1 / x) Exp[-1/x]
R\[ScriptCapitalM]f[y_] := (1 - y) Exp[1/y] ExpIntegralEi[-1/y] / (π y)
R\[ScriptCapitalM][2][B_, θc_, M0_][θ_] := - M0 (R\[ScriptCapitalM]f[B(θc - θ)] - R\[ScriptCapitalM]f[B(θc + θ)])
ruleB[2][f_, g_] := π / (2 1.3578383417066) (Δ g[θc] f[θc]^(-Δ[2]-1) f'[θc] - g'[θc] / f[θc]^Δ[2])
eqLow[D_:2][f_, g_][m_] := SeriesCoefficient[
R\[ScriptCapitalM][D][B, θc, M0][θ] + f[θ]^β[D] Gl'[g[θ] / f[θ]^Δ[D]],
{θ, θc, m},
Assumptions -> Join[$Assumptions, {θ < θc, θ > 1}]
]
eqHigh[D_:2][f_, g_][m_] := SeriesCoefficient[
R\[ScriptCapitalM][D][B, θc, M0][θ] + (-f[θ])^β[D] Gh'[g[θ] / (-f[θ])^Δ[D]],
{θ, 0, m},
Assumptions -> Join[$Assumptions, {θ > 0, θ < 1}]
]
eqMid[D_:2][f_, g_][m_] := SeriesCoefficient[
R\[ScriptCapitalM][D][B, θc, M0][θ] + g[θ]^(1/δ[D]) ((2 - α[D]) Φ'[η] - η Φ'[η]) / Δ[D] /. η -> f[θ] / g[θ]^(1 / Δ[D]),
{θ, 1, m},
Assumptions -> Join[$Assumptions, {θ > 0, θ < θc}]
]
Φs = {
-1.197733383797993,
-0.318810124891,
0.110886196683,
0.01642689465,
-2.639978 10^-4,
-5.140526 10^-4,
2.08856 10^-4,
-4.4819 10^-5
}
Gls = {
0,
-1.3578383417066,
-0.048953289720,
0.038863932,
-0.068362119,
0.18388370,
-0.6591714,
2.937665,
-15.61
}
Ghs = {
0,
0,
-1.8452280782328,
0,
8.333711750,
0,
-95.16896,
0,
1457.62,
0,
-25891
}
dRule[sym_][f_, i_] := Derivative[i[[1]] - 1][sym][0] -> f (i[[1]] - 1)!
ΦRules = MapIndexed[dRule[Φ], Φs];
GlRules = MapIndexed[dRule[Gl], Gls];
GhRules = MapIndexed[dRule[Gh], Ghs];
rules[f_, g_] := Join[{B -> ruleB[2][f, g]}, ΦRules, GlRules, GhRules]
eq[D_:2][f_, g_][m_] := Select[Flatten[{eqLow[D][f, g][#], eqMid[D][f, g][#], eqHigh[D][f, g][#]} & /@ Range[0, m]] /. rules[f, g], # != 0 &]
EndPackage[]
|
