Markiyan Hirnyk

Markiyan Hirnyk
10 years, 161 days


These are answers submitted by Markiyan Hirnyk

Partial answer

April 21 2016 Markiyan Hirnyk 6899

can be done as follows. We find the rational parametrization in t  of the parabola and note that it is enough t=9*n with an arbitrary integer n to this end.

 

with(algcurves):````

L := parametrization(4*x^2-12*x*y+9*y^2-7*x+9*y+369, x, y, t)

[(1/9)*t^2-(1/3)*t+367, (2/27)*t^2-(1/9)*t+244]

(1)

x := L[1]

(1/9)*t^2-(1/3)*t+367

(2)

y := L[2]

(2/27)*t^2-(1/9)*t+244

(3)

eval([x, y], t = 9*n)

[9*n^2-3*n+367, 6*n^2-n+244]

(4)

``

``

Unfortunately, the command

>isolve({(1/9)*t^2-(1/3)*t+367 = n, (2/27)*t^2-(1/9)*t+244 = k});

produces an incorrect answer

{k = (2/19683)*(4779+(1/1162261467)*_Z1)^2+673/3-(1/282429536481)*_Z1, n = (1/6561)*(4779+(1/1162261467)*_Z1)^2+308-(1/94143178827)*_Z1, t = 177+(1/31381059609)*_Z1}

Mathematica finds the complete answer

(C[1] \[Element] Integers && x == 369 - 9 C[1] + 9 C[1]^2 &&
y == 246 - 7 C[1] + 6 C[1]^2) || (C[1] \[Element] Integers &&
x == 367 + 3 C[1] + 9 C[1]^2 && y == 244 + C[1] + 6 C[1]^2)

 

Download integer_points_on_parabola.mw

General case

April 12 2016 Markiyan Hirnyk 6899

Following Igor Hlivka, the two-dimensional normal distribution is created and its characteristics are found. Specifying the parameters, a sample of the size 1000 is created and plotted. See section

 Drawing values from the distribution

in the Wiki article for math background.



restart; with(Statistics); with(LinearAlgebra); with(plots)

h := proc (i) options operator, arrow; x[i] end proc; g := proc (i) options operator, arrow; mu[i] end proc; w := proc (i) options operator, arrow; sigma[i] end proc; N := 2; X := Vector(N, h); M := Vector(N, g); S := Vector(N, w); CorrMat := Matrix([seq([seq(`if`(i = j, 1, rho), j = 1 .. N)], i = 1 .. N)]); C1 := simplify(S.Transpose(S)); CoVar := Zip(`*`, C1, CorrMat)

X := Vector(2, {(1) = x[1], (2) = x[2]})

 

M := Vector(2, {(1) = mu[1], (2) = mu[2]})

 

S := Vector(2, {(1) = sigma[1], (2) = sigma[2]})

 

CorrMat := Matrix(2, 2, {(1, 1) = 1, (1, 2) = rho, (2, 1) = rho, (2, 2) = 1})

 

CoVar := Matrix(2, 2, {(1, 1) = sigma[1]^2, (1, 2) = sigma[1]*sigma[2]*rho, (2, 1) = sigma[1]*sigma[2]*rho, (2, 2) = sigma[2]^2})

(1)

p1 := exp(-(1/2)*Typesetting:-delayDotProduct(Typesetting:-delayDotProduct(LinearAlgebra:-Transpose(X-M), 1/CoVar), X-M))/((2*Pi)^((1/2)*N)*sqrt(Determinant(CoVar))):

bipdf := simplify(p1, symbolic);

(1/2)*exp((1/2)*((x[2]-mu[2])^2*sigma[1]^2-2*rho*sigma[2]*(x[2]-mu[2])*(x[1]-mu[1])*sigma[1]+sigma[2]^2*(x[1]-mu[1])^2)/(sigma[1]^2*sigma[2]^2*(rho-1)*(rho+1)))/(sigma[1]*sigma[2]*(-rho^2+1)^(1/2)*Pi)

(2)

MarDist[X[1]] := `assuming`([int(bipdf, x[2] = -infinity .. infinity)], [-1 < rho and rho < 1, sigma[1] > 0, sigma[2] > 0]);

(1/2)*2^(1/2)*exp(-(1/2)*(mu[1]^2-2*mu[1]*x[1]+x[1]^2)/sigma[1]^2)/(Pi^(1/2)*sigma[1])

(3)

MarDist[X[2]] := `assuming`([int(bipdf, x[1] = -infinity .. infinity)], [-1 < rho and rho < 1, sigma[1] > 0, sigma[2] > 0]);

(1/2)*2^(1/2)*exp(-(1/2)*(mu[2]^2-2*mu[2]*x[2]+x[2]^2)/sigma[2]^2)/(sigma[2]*Pi^(1/2))

(4)

Conditional PDF of  "X|Y = y"          "f[X|Y]= (f[X,Y](x,y))/(f[Y](y))"
CondDist[XY] := simplify(bipdf/MarDist[X[2]], symbolic)

(1/2)*2^(1/2)*exp((1/2)*(sigma[1]*(x[2]-mu[2])*rho-sigma[2]*(x[1]-3))^2/(sigma[1]^2*sigma[2]^2*(rho-1)*(rho+1)))/((-rho^2+1)^(1/2)*Pi^(1/2)*sigma[1])

(5)

mu[1] := 3:``

mu[2] := 15: 

sigma[1] := 1:

rho := -.9:NULL

CoVar

Matrix([[sigma[1]^2, sigma[1]*sigma[2]*rho], [sigma[1]*sigma[2]*rho, sigma[2]^2]])

(6)

A := LUDecomposition(CoVar, method = Cholesky)

Matrix([[1, 0], [-1.8, .8717797887]])

(7)

``

MU := Vector[row]([mu[1], mu[2]])

Vector[row]([3, 15])

(8)

sx := Sample(Normal(0, 1), 1000)

Vector(4, {(1) = ` 1 .. 1000 `*Vector[row], (2) = `Data Type: `*float[8], (3) = `Storage: `*rectangular, (4) = `Order: `*Fortran_order})

(9)

sy := Sample(Normal(0, 1), 1000)

Vector(4, {(1) = ` 1 .. 1000 `*Vector[row], (2) = `Data Type: `*float[8], (3) = `Storage: `*rectangular, (4) = `Order: `*Fortran_order})

(10)

sample := zip(proc (x, y) options operator, arrow; [x, y] end proc, sx, sy)

Vector(4, {(1) = ` 1 .. 1000 `*Vector[row], (2) = `Data Type: `*anything, (3) = `Storage: `*rectangular, (4) = `Order: `*Fortran_order})

(11)

B := map(proc (c) options operator, arrow; MU+convert(Multiply(A, convert(c, Vector[column])), Vector[row]) end proc, sample)

Vector(4, {(1) = ` 1 .. 1000 `*Vector[row], (2) = `Data Type: `*anything, (3) = `Storage: `*rectangular, (4) = `Order: `*Fortran_order})

(12)
pointplot(map(proc (c) options operator, arrow; convert(c, list) end proc, B), scaling = constrained);

 

pointplot(map(proc (c) options operator, arrow; convert(c, list) end proc, B), scaling = constrained)

 

``



Download _adjusted_sample_1000.mw

 

 

 

 



 

It appears there is a lot of such polynomials :

 

bc := {((D@@1)(g))(1)+(1/2)*((D@@1)(s))(1)^2 = 0, g(0) = 0, s(0) = 0, ((D@@1)(s))(0) = 0, ((D@@2)(s))(1) = 0, ((D@@3)(s))(1) = 0}

{(D(g))(1)+(1/2)*(D(s))(1)^2 = 0, g(0) = 0, s(0) = 0, (D(s))(0) = 0, ((D@@2)(s))(1) = 0, ((D@@3)(s))(1) = 0}

(1)

dsys := {diff(g(x), `$`(x, 4)) = 0, diff(s(x), `$`(x, 5)) = 0}

{diff(diff(diff(diff(diff(s(x), x), x), x), x), x) = 0, diff(diff(diff(diff(g(x), x), x), x), x) = 0}

(2)

sol := dsolve(`union`(dsys, bc), {g(x), s(x)})

{g(x) = (1/6)*(-_C1^2-2*_C3-2*_C4)*x^3+(1/2)*_C3*x^2+_C4*x, s(x) = (1/4)*_C1*x^4-_C1*x^3+(3/2)*_C1*x^2}

(3)

eval(sol, [_C1 = 1, _C3 = 1, _C4 = -1])

{g(x) = -(1/6)*x^3+(1/2)*x^2-x, s(x) = (1/4)*x^4-x^3+(3/2)*x^2}

(4)

``

 

Download polynomials.mw

 

Another way

April 06 2016 Markiyan Hirnyk 6899

 

match((a+b+1)*t = -a+b+1, t, 's')

true

(1)

s

{a = 0, b = -1}

(2)

``


Look in ?match for info.

Download another_way.mw

PS. The match command works in several dimensions too, whereas the solve,identity command works in one-dimensional case only.

Specifying the parameters, it can be solved with the DirectSearch package which should be downloaded from http://www.maplesoft.com/applications/view.aspx?SID=101333 and installed in your Maple >=13.

 

restart; A := eval(-(((beta*eta^2-(1/2)*beta+1)*p^2-(1/2)*beta^3)*KummerM(((-beta+2)*p^2-beta^3)/(4*p^2), 1, beta*eta^2)+(1/2)*KummerM(((-beta+6)*p^2-beta^3)/(4*p^2), 1, beta*eta^2)*((beta-2)*p^2+beta^3))*exp(-(1/2)*beta*eta^2)/(p^2*beta) = 0, [p = 1, eta = 1]);

-(((1/2)*beta+1-(1/2)*beta^3)*KummerM(-(1/4)*beta^3-(1/4)*beta+1/2, 1, beta)+(1/2)*KummerM(-(1/4)*beta^3-(1/4)*beta+3/2, 1, beta)*(beta^3+beta-2))*exp(-(1/2)*beta)/beta = 0

(1)

DirectSearch:-SolveEquations([Re(eval(A, beta = a+I*b)), Im(eval(A, beta = a+I*b))], {a = -5 .. 5, b = -5 .. 5}, AllSolutions, number = 300,solutions=10);

Matrix(10, 4, {(1, 1) = 0.1948977161e-26, (1, 2) = Vector(2, {(1) = HFloat(-1.0451202215029108e-14), (2) = HFloat(-4.289230156047799e-14)}), (1, 3) = [a = -4.39157583457543, b = 0.6122407479e-14], (1, 4) = 80, (2, 1) = 0.2468683635e-26, (2, 2) = Vector(2, {(1) = HFloat(4.76220423508406e-14), (2) = HFloat(-1.4171263788954476e-14)}), (2, 3) = [a = 0.9353300437e-13, b = .705280177757183], (2, 4) = 80, (3, 1) = 0.1107221646e-25, (3, 2) = Vector(2, {(1) = HFloat(1.0344270982184144e-13), (2) = HFloat(1.928269298431157e-14)}), (3, 3) = [a = 0.2043030414e-12, b = -.705280177757170], (3, 4) = 112, (4, 1) = 0.1498604938e-25, (4, 2) = Vector(2, {(1) = HFloat(-2.111504886120853e-14), (2) = HFloat(-1.2058276865205994e-13)}), (4, 3) = [a = 4.01763280323693, b = -0.1881620268e-13], (4, 4) = 82, (5, 1) = 0.2098144474e-25, (5, 2) = Vector(2, {(1) = HFloat(1.3962856649975636e-13), (2) = HFloat(-3.853969581723759e-14)}), (5, 3) = [a = -0.2730763486e-13, b = -3.24247382580224], (5, 4) = 87, (6, 1) = 0.5674455128e-25, (6, 2) = Vector(2, {(1) = HFloat(-7.700514667200426e-14), (2) = HFloat(2.2542129150220777e-13)}), (6, 3) = [a = 1.87588952396935, b = 0.7554868781e-13], (6, 4) = 81, (7, 1) = 0.8186683329e-25, (7, 2) = Vector(2, {(1) = HFloat(2.5763416746195455e-13), (2) = HFloat(-1.244647301461954e-13)}), (7, 3) = [a = 0.5969031808e-13, b = 2.71279028626207], (7, 4) = 54, (8, 1) = 0.9325892563e-25, (8, 2) = Vector(2, {(1) = HFloat(5.836258384818908e-14), (2) = HFloat(2.9975445690543977e-13)}), (8, 3) = [a = -1.87588952396936, b = 0.1004610334e-12], (8, 4) = 80, (9, 1) = 0.3667927203e-24, (9, 2) = Vector(2, {(1) = HFloat(-1.582465938541118e-13), (2) = HFloat(-5.845945054753967e-13)}), (9, 3) = [a = 0.4867515661e-13, b = -2.04633245164591], (9, 4) = 74, (10, 1) = 0.3928454083e-24, (10, 2) = Vector(2, {(1) = HFloat(-2.3898361217675557e-13), (2) = HFloat(5.794240600653144e-13)}), (10, 3) = [a = 0.4650368771e-13, b = 3.24247382580213], (10, 4) = 89})

(2)

``

Ten roots out of 29 roots in the square are found. In order to get more solutions, one should increase the value of the solutions option.

Download roots.mw

One should exploit the Newton-Leibnitz formula:


J := `assuming`([int(t^n*ln(t)^n, t = 0 .. x)], [n::posint, x > 0]);

-(n^2*GAMMA(n, -(n+1)*ln(x))*(-n-1)^(-n-1)+n*GAMMA(n, -(n+1)*ln(x))*(-n-1)^(-n-1)-x*(ln(x)*x)^n)/(n+1)

(1)

`assuming`([simplify(diff(J, x)-x^n*ln(x)^n)], [n::posint, x > 0])

0

(2)

``


Download direct_calculation.mw

Improvement

April 02 2016 Markiyan Hirnyk 6899

In fact, this is an improvement of the sum command: see ?sum for info. There is  _EnvFormal option to this end. For example,

 

restart; S := convert(1/(3*x+2), FPS, x)

Sum((1/2)*2^(-k)*(-3)^k*x^k, k = 0 .. infinity)

(1)

value(S)

sum((1/2)*2^(-k)*(-3)^k*x^k, k = 0 .. infinity)

(2)

_EnvFormal := true

true

(3)

value(S)

1/(3*x+2)

(4)

``

 

Download improvement.mw

Another way is 

>restart;sum((1/2)*2^(-k)*(-3)^k*x^k, k = 0 .. infinity)assuming abs(x)<1/ 6;

1/(3*x+2)

The DirectSearch package should be downloaded from http://www.maplesoft.com/applications/view.aspx?SID=101333 and installed in your Maple (>=13). Here is my solution done with the DirectSearch. A more powerful comp is required to obtain more points by increasing the number option.

 

restart; f1 := cos(t1)+1.35*cos(p1)-x; f2 := sin(t1)+1.35*sin(p1)-y; f3 := cos(t2)+1.35*cos(p2)-x+1.15; f4 := sin(t2)+1.35*sin(p2)-y; f5 := cos(t3)+1.35*cos(p3)-x+.575; f6 := sin(t3)+1.35*sin(p3)-y+.995; f7 := cos(t1)*cos(t2)*cos(t3)*sin(p1)*sin(p2)*sin(p3)-cos(t1)*cos(t2)*sin(t3)*sin(p1)*sin(p2)*cos(p3)-cos(t1)*sin(t2)*cos(t3)*sin(p1)*cos(p2)*sin(p3)+cos(t1)*sin(t2)*sin(t3)*sin(p1)*cos(p2)*cos(p3)-sin(t1)*cos(t2)*cos(t3)*cos(p1)*sin(p2)*sin(p3)+sin(t1)*cos(t2)*sin(t3)*cos(p1)*sin(p2)*cos(p3)+sin(t1)*sin(t2)*cos(t3)*cos(p1)*cos(p2)*sin(p3)-sin(t1)*sin(t2)*sin(t3)*cos(p1)*cos(p2)*cos(p3); f8 := cos(t1)*sin(p1)-sin(t1)*cos(p1) >= 0; f9 := cos(t2)*sin(p2)-sin(t2)*cos(p2) >= 0; f10 := cos(t3)*sin(p3)-sin(t3)*cos(p3) >= 0

sol := DirectSearch:-SolveEquations([seq(f || j, j = 1 .. 7)], {f10, f8, f9, p1 = -Pi .. Pi, p2 = -Pi .. Pi, p3 = -Pi .. Pi, t1 = -Pi .. Pi, t2 = -Pi .. Pi, t3 = -Pi .. Pi, x = -10 .. 10, y = -10 .. 10}, AllSolutions, method = quadratic, number = 200);

RTABLE(18446744074192502782, anything, Matrix, rectangular, Fortran_order, [], 2, 1 .. 182, 1 .. 4)

(1)

sol[1];

Vector[row](4, {(1) = 0.1071032753e-27, (2) = Vector(7, {(1) = HFloat(-2.55351295663786e-15), (2) = HFloat(6.661338147750939e-16), (3) = HFloat(4.884981308350689e-15), (4) = HFloat(-3.1086244689504383e-15), (5) = HFloat(7.105427357601002e-15), (6) = HFloat(3.9968028886505635e-15), (7) = HFloat(3.885780586188048e-16)}), (3) = [p1 = 2.81233098296781, p2 = 2.54814516295493, p3 = -2.62961132730027, t1 = 1.07095279467161, t2 = 2.54814516295493, t3 = 1.76837387160141, x = -.798191584527774, y = 1.31417257237213], (4) = 1032})

(2)

``

 

 

plots:-pointplot3d([seq([rhs(sol(1 .. 53, 3)[j][7]), rhs(sol(1 .. 53, 3)[j][8]), rhs(sol(1 .. 53, 3)[j][4])], j = 1 .. 53)], axes = frame)

 

``

 

Download By_DS.mw

 

 

Link

March 15 2016 Markiyan Hirnyk 6899

Hope this will be useful for you.

This can be done as follows.


series(t^2*exp(-t)/ln(t), t = 1, 3)

series(exp(-1)/(t-1)+(3/2)*exp(-1)+O((t-1)^1),t = 1,1)

(1)

A := int(t^2*exp(-t)/ln(t), t = 0 .. 1/2, numeric)

-0.2953824719e-1

(2)

B := int(t^2*exp(-t)/ln(t), t = 3/2 .. infinity, numeric)

1.645313668

(3)

C := int(t^2*exp(-t)/ln(t)-exp(-1)/(t-1), t = 1/2 .. 3/2, numeric)

.5377722704

(4)

E := int(exp(-1)/(t-1), t = 1/2 .. 3/2, CPV)

0

(5)

A+B+C+E

2.153547691

(6)

``


Download 5steps.mw

Somewhat extending the Axel's suggestion, one can use

 

eqs := {-a*y+x, -a*x+y};

piecewise(x <> 0, piecewise(a-1 <> 0, [], a+1 <> 0, [], a = 0, [], And(a-1 = 0, a+1 = 0, a <> 0), [{y = x/a}]), a-1 <> 0, piecewise(x <> 0, [], x = 0, [{y = x/a}]), a+1 <> 0, piecewise(x <> 0, [], x = 0, [{y = x/a}]), a = 0, piecewise(x = 0, [{y = 0}], x <> 0, []), And(x = 0, a-1 = 0, a+1 = 0, a <> 0), [{y = x/a}])

(1)

``

 

Download parametric=full.mw

Too

March 05 2016 Markiyan Hirnyk 6899

J := int(arctan(x)*ln(1+1/x^2), x = 0 .. infinity);
IntegrationTools:-Parts(J, ln(1+1/x^2));

                          (1/6)*Pi^2 

Another way

March 02 2016 Markiyan Hirnyk 6899

This can be directly done by

 

`assuming`([solve(abs(n^2/(n^2+31*n+228)-1) < epsilon, n)], [epsilon > 0])

{n < -(1/2)*(31*epsilon+31+(49*epsilon^2+1010*epsilon+961)^(1/2))/epsilon}, {n < -228/31, (1/2)*(-31*epsilon-31+(49*epsilon^2+1010*epsilon+961)^(1/2))/epsilon < n}

(1)

``


One can take

N=

Download N.mw

Reference

March 02 2016 Markiyan Hirnyk 6899

Look in http://www.maplesoft.com/applications/view.aspx?SID=33406

and http://www.mapleprimes.com/posts/89226-WithOrthogonalExpansions

These links can be found by the "fourier" search in MaplePrimes at the top of this page.

 

1 2 3 4 5 6 7 Last Page 1 of 134