Question: How to add vector fields

Hallo every body 

How to add vector fields to the figure of this example of a three-dimensional differential system.

in maple 18

Porgram_of_corollary_1_in_Maple.mw

NULL

restart

X[j] := x^3*a[j, 0]+x^2*y*a[j, 1]+x^2*z*a[j, 2]+x*y^2*a[j, 3]+x*y*z*a[j, 4]+x*z^2*a[j, 5]+y^3*a[j, 6]+y^2*z*a[j, 7]+y*z^2*a[j, 8]+z^3*a[j, 9]

x^3*a[j, 0]+x^2*y*a[j, 1]+x^2*z*a[j, 2]+x*y^2*a[j, 3]+x*y*z*a[j, 4]+x*z^2*a[j, 5]+y^3*a[j, 6]+y^2*z*a[j, 7]+y*z^2*a[j, 8]+z^3*a[j, 9]

(1)

s := sum(epsilon^j*X[j], j = 0 .. 2)

x^3*a[0, 0]+x^2*y*a[0, 1]+x^2*z*a[0, 2]+x*y^2*a[0, 3]+x*y*z*a[0, 4]+x*z^2*a[0, 5]+y^3*a[0, 6]+y^2*z*a[0, 7]+y*z^2*a[0, 8]+z^3*a[0, 9]+epsilon*(x^3*a[1, 0]+x^2*y*a[1, 1]+x^2*z*a[1, 2]+x*y^2*a[1, 3]+x*y*z*a[1, 4]+x*z^2*a[1, 5]+y^3*a[1, 6]+y^2*z*a[1, 7]+y*z^2*a[1, 8]+z^3*a[1, 9])+epsilon^2*(x^3*a[2, 0]+x^2*y*a[2, 1]+x^2*z*a[2, 2]+x*y^2*a[2, 3]+x*y*z*a[2, 4]+x*z^2*a[2, 5]+y^3*a[2, 6]+y^2*z*a[2, 7]+y*z^2*a[2, 8]+z^3*a[2, 9])

(2)

s1 := subs(a = b, s)

x^3*b[0, 0]+x^2*y*b[0, 1]+x^2*z*b[0, 2]+x*y^2*b[0, 3]+x*y*z*b[0, 4]+x*z^2*b[0, 5]+y^3*b[0, 6]+y^2*z*b[0, 7]+y*z^2*b[0, 8]+z^3*b[0, 9]+epsilon*(x^3*b[1, 0]+x^2*y*b[1, 1]+x^2*z*b[1, 2]+x*y^2*b[1, 3]+x*y*z*b[1, 4]+x*z^2*b[1, 5]+y^3*b[1, 6]+y^2*z*b[1, 7]+y*z^2*b[1, 8]+z^3*b[1, 9])+epsilon^2*(x^3*b[2, 0]+x^2*y*b[2, 1]+x^2*z*b[2, 2]+x*y^2*b[2, 3]+x*y*z*b[2, 4]+x*z^2*b[2, 5]+y^3*b[2, 6]+y^2*z*b[2, 7]+y*z^2*b[2, 8]+z^3*b[2, 9])

(3)

s2 := subs(a = c, s)

x^3*c[0, 0]+x^2*y*c[0, 1]+x^2*z*c[0, 2]+x*y^2*c[0, 3]+x*y*z*c[0, 4]+x*z^2*c[0, 5]+y^3*c[0, 6]+y^2*z*c[0, 7]+y*z^2*c[0, 8]+z^3*c[0, 9]+epsilon*(x^3*c[1, 0]+x^2*y*c[1, 1]+x^2*z*c[1, 2]+x*y^2*c[1, 3]+x*y*z*c[1, 4]+x*z^2*c[1, 5]+y^3*c[1, 6]+y^2*z*c[1, 7]+y*z^2*c[1, 8]+z^3*c[1, 9])+epsilon^2*(x^3*c[2, 0]+x^2*y*c[2, 1]+x^2*z*c[2, 2]+x*y^2*c[2, 3]+x*y*z*c[2, 4]+x*z^2*c[2, 5]+y^3*c[2, 6]+y^2*z*c[2, 7]+y*z^2*c[2, 8]+z^3*c[2, 9])

(4)

Considérons le système suivant:

a[1] := 0; c[1] := 0

a[0, 9] := 0; c[0, 8] := 0; b[0, 7] := 0; a[0, 4] := 0; a[0, 7] := 0; c[0, 3] := 0; c[0, 0] := 0; c[0, 5] := 0; b[0, 4] := 0; a[0, 2] := 0; c[0, 6] := 0; c[0, 1] := 0; c[0, 7] := 0; a[0, 8] := 0; b[0, 5] := 0

b0 := 5; a[4] := 0; c[4] := 0; c[2, 9] := 0; c[2, 2] := 0; c[2, 7] := 0; a[2, 5] := 0; b[2, 8] := 0; a[2, 0] := 0; b[2, 6] := 0; b[2, 1] := 0; a[2, 3] := 0; b[0, 9] := 0

b[1] := 0; b[2] := 0; b[3] := 0; b[4] := 0; a[1, 2] := 0; a[1, 1] := 0; a[1, 4] := 0; a[1, 6] := 0; a[1, 7] := 0; a[1, 8] := 0; a[1, 9] := 0; a[2, 9] := 0; a[2, 8] := 0; a[2, 7] := 0; a[2, 6] := 0; a[2, 4] := 0; a[2, 2] := 0; a[2, 1] := 0; b[1, 0] := 0; b[1, 2] := 0; b[1, 3] := 0; b[1, 4] := 0; b[1, 5] := 0; b[1, 7] := 0; b[1, 9] := 0; b[2, 0] := 0; b[2, 2] := 0; b[2, 3] := 0; b[2, 4] := 0; b[2, 5] := 0; b[2, 7] := 0; b[2, 9] := 0; c[1, 0] := 0; c[1, 1] := 0; c[1, 3] := 0; c[1, 4] := 0; c[1, 5] := 0; c[1, 6] := 0; c[1, 8] := 0; c[2, 0] := 0; c[2, 1] := 0; c[2, 3] := 0; c[2, 4] := 0; c[2, 5] := 0; c[2, 6] := 0; c[2, 8] := 0; b[0, 2] := 0; c[1, 7] := 0

a[1, 0] := 0; a[1, 3] := 0; a[1, 5] := 0; b[1, 1] := 0; b[1, 6] := 0; b[1, 8] := 0; c[1, 2] := 0; c[1, 9] := 0; a[3] := 0; c[3] := 0; a[2] := 1/2; c[2] := 3/2; a[0, 0] := -1/2; a[0, 3] := 5/4; a[0, 1] := 0; a[0, 5] := 0; a[0, 6] := 0; b[0, 6] := -1; b[0, 1] := 3/2; b[0, 0] := 0; b[0, 3] := 0; b[0, 8] := 0; c[0, 2] := 0; c[0, 4] := 0; c[0, 9] := -1/3

eq1 := (epsilon^4*a[4]+epsilon^3*a[3]+epsilon^2*a[2]+epsilon*a[1])*x-(epsilon^4*b[4]+epsilon^3*b[3]+epsilon^2*b[2]+epsilon*b[1]+b0)*y+s

(1/2)*epsilon^2*x-5*y-(1/2)*x^3+(5/4)*x*y^2

(5)

eq2 := (epsilon^4*b[4]+epsilon^3*b[3]+epsilon^2*b[2]+epsilon*b[1]+b0)*x+(epsilon^4*a[4]+epsilon^3*a[3]+epsilon^2*a[2]+epsilon*a[1])*y+s1

5*x+(1/2)*epsilon^2*y+(3/2)*x^2*y-y^3

(6)

eq3 := (epsilon^4*c[4]+epsilon^3*c[3]+epsilon^2*c[2]+epsilon*c[1])*z+s2

(3/2)*epsilon^2*z-(1/3)*z^3

(7)

Faisons le changement (x,y,z)=(εX,εY,εZ)

 

x := epsilon*X; y := epsilon*Y; z := epsilon*Z

epsilon*X

 

epsilon*Y

 

epsilon*Z

(8)

Xpoint := collect(eq1/epsilon, epsilon)

((1/2)*X-(1/2)*X^3+(5/4)*X*Y^2)*epsilon^2-5*Y

(9)

Ypoint := collect(eq2/epsilon, epsilon)

((1/2)*Y+(3/2)*X^2*Y-Y^3)*epsilon^2+5*X

(10)

Zpoint := collect(eq3/epsilon, epsilon)

((3/2)*Z-(1/3)*Z^3)*epsilon^2

(11)

Faisons le changement (X, Y, Z) = (`ϱ`*cos(theta), `ϱ`*sin(theta), eta)

 

X := `ϱ`*cos(theta); Y := `ϱ`*sin(theta); Z := eta

`ϱ`*cos(theta)

 

`ϱ`*sin(theta)

 

eta

(12)

`ϱt` := collect(simplify((X*Xpoint+Y*Ypoint)/`ϱ`), epsilon)

-(1/4)*`ϱ`*epsilon^2*(17*`ϱ`^2*cos(theta)^4-19*cos(theta)^2*`ϱ`^2+4*`ϱ`^2-2)

(13)

`θt` := collect(simplify((X*Ypoint-Xpoint*Y)/`ϱ`^2), epsilon)

5+((17/4)*`ϱ`^2*cos(theta)^3*sin(theta)-(9/4)*`ϱ`^2*sin(theta)*cos(theta))*epsilon^2

(14)

`ηt` := collect(Zpoint, epsilon)

((3/2)*eta-(1/3)*eta^3)*epsilon^2

(15)

Utilisons le développpement de taylor

p := series(`ϱt`/`θt`, epsilon, 5)

series(-((1/20)*`ϱ`*(17*`ϱ`^2*cos(theta)^4-19*cos(theta)^2*`ϱ`^2+4*`ϱ`^2-2))*epsilon^2+((1/100)*`ϱ`*(17*`ϱ`^2*cos(theta)^4-19*cos(theta)^2*`ϱ`^2+4*`ϱ`^2-2)*((17/4)*`ϱ`^2*cos(theta)^3*sin(theta)-(9/4)*`ϱ`^2*sin(theta)*cos(theta)))*epsilon^4+O(epsilon^6),epsilon,6)

(16)

q := series(`ηt`/`θt`, epsilon, 5)

series(((3/10)*eta-(1/15)*eta^3)*epsilon^2+((1/5)*(-(3/10)*eta+(1/15)*eta^3)*((17/4)*`ϱ`^2*cos(theta)^3*sin(theta)-(9/4)*`ϱ`^2*sin(theta)*cos(theta)))*epsilon^4+O(epsilon^6),epsilon,6)

(17)

NULL

Averaging d'ordre 1

Les fonctions F11 et F21 sont données comme suit:

NULL

F11 := coeff(p, epsilon)

0

(18)

F21 := coeff(q, epsilon)

0

(19)

NULL

Calculons les fonctions moyennées f11et f12

f11 := (int(F11, theta = 0 .. 2*Pi))/(2*Pi)

0

(20)

f12 := (int(F21, theta = 0 .. 2*Pi))/(2*Pi)

0

(21)

solve({f11 = 0, f12 = 0}, {eta, `ϱ`})

{eta = eta, `ϱ` = `ϱ`}

(22)

NULL

Averaging d'ordre 2

NULL

F12 := simplify(coeff(p, epsilon^2))

-(1/20)*`ϱ`*(17*`ϱ`^2*cos(theta)^4-19*cos(theta)^2*`ϱ`^2+4*`ϱ`^2-2)

(23)

F22 := simplify(coeff(q, epsilon^2))

(3/10)*eta-(1/15)*eta^3

(24)

NULL

Calculons les fonctions moyennées "f21 "et "f22"

f21 := simplify((int(F12, theta = 0 .. 2*Pi))/(2*Pi))

-(1/160)*`ϱ`*(7*`ϱ`^2-16)

(25)

f22 := simplify((int(F22, theta = 0 .. 2*Pi))/(2*Pi))

-(1/30)*eta*(2*eta^2-9)

(26)

solve({f21 = 0, f22 = 0}, {eta, `ϱ`})

{eta = 0, `ϱ` = 0}, {eta = 3*RootOf(2*_Z^2-1), `ϱ` = 0}, {eta = 0, `ϱ` = 4*RootOf(7*_Z^2-1)}, {eta = 3*RootOf(2*_Z^2-1), `ϱ` = 4*RootOf(7*_Z^2-1)}

(27)

allvalues({eta = 0, `ϱ` = 4*RootOf(7*_Z^2-1)})

{eta = 0, `ϱ` = (4/7)*7^(1/2)}, {eta = 0, `ϱ` = -(4/7)*7^(1/2)}

(28)

allvalues({eta = 3*RootOf(2*_Z^2-1), `ϱ` = 4*RootOf(7*_Z^2-1)})

{eta = (3/2)*2^(1/2), `ϱ` = (4/7)*7^(1/2)}, {eta = -(3/2)*2^(1/2), `ϱ` = (4/7)*7^(1/2)}, {eta = (3/2)*2^(1/2), `ϱ` = -(4/7)*7^(1/2)}, {eta = -(3/2)*2^(1/2), `ϱ` = -(4/7)*7^(1/2)}

(29)

NULL

with(VectorCalculus)

M, d := Jacobian([f21, f22], [`ϱ`, eta] = [(4/7)*sqrt(7), 0], 'determinant')

Matrix(%id = 18446744074358842782), -3/50

(30)

factor(d)

-3/50

(31)

M1, d1 := Jacobian([f21, f22], [`ϱ`, eta] = [(4/7)*sqrt(7), (3/2)*sqrt(2)], 'determinant')

Matrix(%id = 18446744074358843142), 3/25

(32)

d1 := factor(d1)

3/25

(33)

M2, d2 := Jacobian([f21, f22], [`ϱ`, eta] = [(4/7)*sqrt(7), -(3/2)*sqrt(2)], 'determinant')

Matrix(%id = 18446744074358843382), 3/25

(34)

factor(d2)

3/25

(35)

restart

with(DEtools):

epsilon := 10^(-2)

1/100

(36)

eq1 := diff(x(t), t) = (1/2)*epsilon^2*x(t)-5*y(t)-(1/2)*x(t)^3+(5/4)*x(t)*y(t)^2

diff(x(t), t) = (1/20000)*x(t)-5*y(t)-(1/2)*x(t)^3+(5/4)*x(t)*y(t)^2

(37)

eq2 := diff(y(t), t) = 5*x(t)+(1/2)*epsilon^2*y(t)+(3/2)*x(t)^2*y(t)-y(t)^3

diff(y(t), t) = 5*x(t)+(1/20000)*y(t)+(3/2)*x(t)^2*y(t)-y(t)^3

(38)

eq3 := diff(z(t), t) = (3/2)*epsilon^2*z(t)-(1/3)*z(t)^3

diff(z(t), t) = (3/20000)*z(t)-(1/3)*z(t)^3

(39)

DEplot3d([eq1, eq2, eq3], [x(t), y(t), z(t)], t = -10 .. 10, [[x(0) = 0.1511857892e-1, y(0) = 0, z(0) = 0], [x(0) = 0.1511857892e-1, y(0) = 0, z(0) = 0.2121320343e-1], [x(0) = 0.1511857892e-1, y(0) = 0, z(0) = -0.2121320343e-1]], linecolor = [blue, red, black], stepsize = 0.1e-1)

 

Download Porgram_of_corollary_1_in_Maple.mw

Please Wait...