## 8662 Reputation

8 years, 348 days

## Functional assignment...

Use functional assignment for the calculating  f   for different points.

See corrected file.

## Implementation of a step-by-step Prim's ...

Unfortunately, GraphTheory:-PrimsAlgorithm command returns a list of edges for the minimal spanning tree not in the same order as these edges are obtained by applying Prim's algorithm. Here is another procedure named  PrimsAlg  that returns a list of edges in the correct order. The procedure does not use any commands from GraphTheory package and works strictly according to the Prim's algorithm.The required parameter  L  is the list of the vertices of the graph given by the coordinates of points on the Euclidean plane. The optional parameter  N   is the number of the vertex with which we begin to build the minimal spanning tree ( by default N=1).

The code of the procedure:

PrimsAlg:=proc(L::listlist, N::posint:=1)
local n, A, dist, L1, V, U, i, p, E;
n:=nops(L);
A:=Matrix(n, (i,j)->evalf(sqrt((L[i,1]-L[j,1])^2+(L[i,2]-L[j,2])^2)));

dist:=proc(U::set(posint), V::set(posint))
local a, u, v, p;
a:=infinity;
for u in U do
for v in V do
if A[u,v]<a then a:=A[u,v]; p:=[u,v,a] fi;
od;
od;
p;
end proc:

L1:=convert(L, set):
V[1]:={N}: U[1]:= {\$1..n} minus V[1];

for i from 1 to n-1 do
p:=dist(V[i],U[i]);
V[i+1]:=V[i] union {p[2]}; U[i+1]:=U[i] minus {p[2]};
E[i]:=p[1..2];
od;

convert(E, list);

end proc:

Example of use for the list of 50 random points in the square  [0,30] x [0,30]:

L:=RandomTools:-Generate(listlist(integer(range = 0 .. 30), 53, 2)):  # Generating 53 random points
M:=convert(convert(L,set)[1..50], list):  #  Selection of 50 different points
M1:=subs({seq(n=M[n], n=1..50)}, PrimsAlg(M)):  # Consecutive list of added edges

M1 := [[[0, 3], [2, 1]], [[0, 3], [2, 7]], [[2, 7], [2, 8]], [[2, 8], [0, 11]], [[2, 8], [4, 11]], [[4, 11], [6, 8]], [[4, 11], [8, 12]], [[8, 12], [9, 16]], [[9, 16], [8, 17]], [[8, 17], [8, 20]], [[8, 20], [8, 21]], [[8, 21], [10, 24]], [[10, 24], [8, 26]], [[8, 26], [8, 29]], [[8, 12], [12, 11]], [[12, 11], [13, 12]], [[8, 17], [4, 16]], [[4, 16], [1, 17]], [[1, 17], [2, 21]], [[8, 26], [4, 25]], [[4, 25], [4, 28]], [[10, 24], [14, 23]], [[14, 23], [14, 22]], [[14, 23], [17, 23]], [[17, 23], [16, 27]], [[16, 27], [15, 28]], [[15, 28], [13, 30]], [[16, 27], [19, 29]], [[17, 23], [21, 22]], [[12, 11], [16, 8]], [[16, 8], [17, 8]], [[17, 8], [18, 8]], [[18, 8], [20, 9]], [[20, 9], [20, 13]], [[20, 13], [22, 13]], [[22, 13], [23, 16]], [[18, 8], [20, 4]], [[20, 4], [18, 2]], [[18, 2], [16, 1]], [[20, 4], [23, 4]], [[23, 4], [27, 4]], [[27, 4], [27, 1]], [[27, 1], [29, 0]], [[27, 1], [24, 0]], [[27, 4], [30, 7]], [[22, 13], [27, 13]], [[23, 16], [28, 19]], [[28, 19], [29, 19]], [[19, 29], [27, 28]]]

Visualization:

A:=plot(M, style=point, color=blue, symbol=solidcircle, symbolsize=15):
B:=seq(plot(M1, color=red, thickness=2), i=1..49):
plots:-display(A, B);
# The minimal spanning tree for M

Addition - animation of the building of the minimal spanning tree  (the code below does not depend on the previous one):

M1 := [[[0, 3], [0, 3]], [[0, 3], [2, 1]], [[0, 3], [2, 7]], [[2, 7], [2, 8]], [[2, 8], [0, 11]], [[2, 8], [4, 11]], [[4, 11], [6, 8]], [[4, 11], [8, 12]], [[8, 12], [9, 16]], [[9, 16], [8, 17]], [[8, 17], [8, 20]], [[8, 20], [8, 21]], [[8, 21], [10, 24]], [[10, 24], [8, 26]], [[8, 26], [8, 29]], [[8, 12], [12, 11]], [[12, 11], [13, 12]], [[8, 17], [4, 16]], [[4, 16], [1, 17]], [[1, 17], [2, 21]], [[8, 26], [4, 25]], [[4, 25], [4, 28]], [[10, 24], [14, 23]], [[14, 23], [14, 22]], [[14, 23], [17, 23]], [[17, 23], [16, 27]], [[16, 27], [15, 28]], [[15, 28], [13, 30]], [[16, 27], [19, 29]], [[17, 23], [21, 22]], [[12, 11], [16, 8]], [[16, 8], [17, 8]], [[17, 8], [18, 8]], [[18, 8], [20, 9]], [[20, 9], [20, 13]], [[20, 13], [22, 13]], [[22, 13], [23, 16]], [[18, 8], [20, 4]], [[20, 4], [18, 2]], [[18, 2], [16, 1]], [[20, 4], [23, 4]], [[23, 4], [27, 4]], [[27, 4], [27, 1]], [[27, 1], [29, 0]], [[27, 1], [24, 0]], [[27, 4], [30, 7]], [[22, 13], [27, 13]], [[23, 16], [28, 19]], [[28, 19], [29, 19]], [[19, 29], [27, 28]]]:
M:={op~(M1)[ ]}:
A:=plot(M, style=point, color=blue, symbol=solidcircle, symbolsize=15):
F:=p->plots:-display(seq(plottools:-line(M1[i][], color=red, thickness=2), i=1..p));
plots:-animate(plots:-display,['F'(p)], p=1..50, background=A, frames=51);

Prim.mw

Edit.

You can use  textplot  command instead of  legend  option to place the legends in any desired place. Of course, this will require a few more efforts.

Example:

with(plots):
A:=plot([sin(x), cos(x)], x=-Pi..3*Pi, color=[red,blue]):
B:=textplot([[3.7,0.8,sin(x)], [-1.8,0.5,cos(x)]], font=[times,bold,
14]):
C:=textplot([[3,1,__, color=red], [-2.5,0.7,__,color=blue]], font=[times,bold,24]):
display(A, B, C, size=[950,200], scaling=constrained);

## PolyhedralSets...

Use  PolyhedralSets  package. See the help on this package for details. Here is an example of the plotting of your first set (all inequalities must be non-strict):

with(PolyhedralSets):
P:=PolyhedralSet({40 <= y, x <= y-40, z <= y-40});
Plot(P, color=khaki, axes=normal, view=[-30..60, 0..90, -40..50], orientation=[40,70]);

We see that this is an unbounded trihedral body.

## animate or Explore...

A lot of syntax errors. Also note that  implicitplot  command only works for two variables. To examine the dependence of the plot on the third variable (I took  p  as the third variable), you can use  animate  or Explore  commands:

restart;
with(plots):
omega[0]:=1: epsilon[1]:=1.04493: epsilon[2]:=0.93259: delta:=0.02:
eq:=(omega[0]*a-a*Omega^2+(3/4)*epsilon[2]*a^3)^2+(delta*Omega*a)^2 = p^2;
animate(implicitplot, [eq, Omega = 0.5 .. 2, a = 0 .. 5, color = blue, thickness = 2, gridrefine = 3], p=0..0.2, frames=120);

# Or
Explore(implicitplot(eq, Omega = 0.5 .. 2, a = 0 .. 5, color = blue, thickness = 2, gridrefine = 3), parameters=[p=0...1.], numframes =120);

Edit.

## Options...

We can noticeably improve the perception if we emphasize the ellipse by its border and/or a color.  fieldstrength=log   option is also useful, which makes it possible to reduce the difference in the size of the arrows:

restart;
with(plots):
f := (x,y) -> `if`(x^2+y^2/2 <= 1, x^2+y^2, 0):
g := (x,y) -> `if`(x^2+y^2/2 <= 1, x+y, 0):
plots:-display(plots:-implicitplot(x^2+y^2/2-1<=0, x=-1.5..1.5,y=-1.5..1.5, thickness=0, filledregions=true, coloring=["LightYellow",white]), fieldplot( [f,g], -1.5..1.5,-1.5..1.5, arrows=SLIM, color="Red", fieldstrength=log));

## Prim's algorithm...

Your algorithm may not provide a minimum spanning tree. You should at each step look for the nearest vertex not for the last vertex of the already constructed tree, but to search for the nearest vertex by scanning all the vertices of the already constructed tree. This will be Prim's algorithm  https://en.wikipedia.org/wiki/Prim%27s_algorithm

Here is an example where your algorithm fails (start from vertex a):

A:=plot([[0,0],[-1,-1],[0,2],[1,3]], style=point, symbol=solidcircle, color=red, symbolsize=20):
B:=plots:-textplot([[0,0,a],[-1,-1,b],[0,2,c],[1,3,d]], color=blue, align={below,right}, font=[times,bold,20]):
plots:-display(A, B, scaling=constrained, axes=none);

Addition. Prim's algorithm is easy to implement in Maple, but now I do not have time to do this (maybe tonight I'll do it)

## timestep and spacestep options...

To improve accuracy, I used timestep and spacestep options. I took the range  t=0..0.5, because already for  t = 0.5  your function is practically equal to 0 (see below). For integration, a spline of 1 order was used :

restart;
k:=5:
EQ:=diff(u(x,t),t)=k*diff(u(x,t), x\$2):
ibc:=u(0,t)=0, u(1,t)=0, u(x,0) = x:
sol:=pdsolve({EQ}, {ibc}, numeric, timestep = 0.001, spacestep = 0.005):
p1:=sol:-plot(u, x=0.5, t=0...0.5, style = line, color = "Blue", legend = "heat Plot", axes=boxed);
M:=op(1, op(1, p1)):
f:=CurveFitting:-Spline(M, t, degree=1):
int(f, t=0..0.5);
# The area under p1
eval(f, t=0.5);  # The value of the function for t=0.5

## extrema, SecondDerivativeTest...

extrema  command by  's'-option  allows you to find all the critical points (in which both partial derivatives are 0). Student[MultivariateCalculus][SecondDerivativeTest]  command allows you to classify these points (LocaMin, LocalMax, or Saddle) :

restart;

f:=(x,y)->x^3 + 3*x*y^2 - 15*x + y^3 - 15*y:
extrema(f(x,y), { }, {x,y}, 's');
# The minimum and maximum values of the function at critical points

s;  # All critical points

for p in s do
p1:=convert(p, list);
print(Student[MultivariateCalculus][SecondDerivativeTest](f(x,y), lhs~(p1)=rhs~(p1)));
od:
# Classification of the critical points

with(plots):
display(plot3d(f(x,y), x=-6..6, y=-6..6, view=-60..60, orientation=[145,60]), pointplot3d(map(t->[rhs~(t)[ ], f(rhs~(t)[ ])]
, s), symbol=solidcircle, color=red, symbolsize=20));  # Plot

Edit.

## In more detail...

f:=x->Statistics:-PDF(Weibull(6.2, 2), x):
A:=plot(f, 0..26, color="DodgerBlue", thickness=3, tickmarks=[spacing(5), spacing(0.02)], axis[2]=[gridlines=8], view=[0..30,0..0.16], size=[600,400], legend=Series1, legendstyle=[font=[times,roman,16], location=right]):
B:=plot([seq([x,f(x)], x=0..26)], style=point, symbol=soliddiamond, symbolsize=18, color="DodgerBlue"):
plots:-display(A, B);

## This is impossible...

Suppose that such a function  F  exists. Then we have for  t=exp(1) :

F(F(exp(t))) = F(t) = F(exp(1))=1  instead of  .

We see that these conditions (F(exp(t)) = t  and  F(F(exp(t))) = 0)

## content, primpart...

y=0.0000125698-0.0000125698*cos(54*x);
a:=convert(%, fraction);
b:=content(rhs(a), cos(54*x));
c:=primpart(rhs(a), cos(54*x));
y=b*``(c);
convert(%, float, 6);

## fsolve...

You just have to solve equation  g(x)=1.2*f(x) :

f:=x->75000*exp(-0.045*x)+100000:
g:=x->300000/(1+100000/3333*exp(-0.08*x)):
fsolve(g(x)=1.2*f(x), x=0..70);
# The result
plot([f, g], 0..70, color=[red,green], legend=[urban,rural]);  # Visualization

## Ways...

Try

```int(evalc(exp(I*x*(2-y))),[x=-infinity..infinity, y=0..2]);  # Or
evalf(Int(evalc(exp(I*x*(2-y))),[x=-infinity..infinity, y=0..2]));
```

## Discrete random variable...

The  F  procedure solves the problem:

F:=(i, x)->rand(i-x .. i+x)( );

Example of use:

seq(F(20, 5), k=1..20);

19, 15, 16, 22, 24, 16, 19, 15, 17, 23, 16, 16, 17, 19, 18, 18, 19, 25, 22, 25

Since there are  i+x-(i-x-1)=2*x+1  numbers in total in the range  i-x .. i+x, the probability of each individual number is  1/(2*x+1) .

Edit.

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