sand15

510 Reputation

10 Badges

4 years, 322 days

MaplePrimes Activity


These are answers submitted by sand15

To complete Carl's answer...

If you are as lazzy as I am, here is maybe a simple trick which could interest you:

Let R a convex region (!!!) :
R can be written by the list L= [P1, ..., PN] of its "corner" points, picked in the clockwise or counter clockwise sense.
For instance a list associated to R1 can be L =  [ [p_l, p_B], [epsilon, p_B], [epsilon,1], [p_l,1] ]

The Trick :

with(PolydrehalSets):
L :=  [ [p_l, p_B], [epsilon, p_B], [epsilon,1], [p_l,1] ]:
obj := PolydrehalSet(L):
rel := Relations(obj)  ; # rel contains the desired relations

PS 1: happily PolyhedraSet works well even on flat polyhedral, that is on polygons

PS2 : unhappily PolyhedraSet needs that each coordinate in L be a rational number
          If it's not so, you need to write L := convert~(L, rational) before calling Polyhedral Set.
          This puts some limitations on this trick, but it can remains usefull...

@Otavio 

From my experience the better way to do an action when clicking a button is to set its attribute onclick this way
onclick = Action( Evaluate( 'function' = 'Myproc', Argument('a1'), ..., Argument('an'))
If necessary an attribute 'target'='T' can be added in the definition of Evaluate(...)
Please see Maplets[Elements](Evaluate] for more details.

Here 'a1', ..., 'an' and 'T' are n elements the "master maplet" must contain (for instance n TextField elements referenced 'a1', ..., 'an' and a Plotter referenced 'T'.

When you clic on the desired button the procedure MyProc will do some operations of arbitrary complexity (just as any Maple's procedure) and run a "slave maplet" named M.

Example: The master maplet contains a TextField element named 'TF' whose field contains the value 10.
You want Myproc to open a maplet with N TextFields (named 'TF1', ..., 'TF10' below)

  1. associate to the opening button the action onclick = Action( Evaluate( 'function' = 'Myproc', Argument('TF'))
     
  2. Define MyProc:
    Myproc := proc(N)
    local M:
    M := Maplet(
        ......
       GridRow( seq( GridCell (TextField[sprintf("TF%d", n)](.....)),  n=1..N) )
       .......
    )
    Maplets[Display](M)
    end proc:
     
  3. I find interesting (a personal opinion), that the master maplet be encapsulated itself in a procedure:
    Main := proc():
    local MyProc, MainMaplet:
    ......

    ... write here the definition of Myproc...
    ... write here the definition of MainMaplet...
    Maplets[Display](MainMaplet);
    ......

    end proc

 

Unfortunately I do not have Maple on this machine: all I wrote hereis the result of a succession of back and forth between this machine and the one I have Maple on... I hope I missed nothing.

Let me know if you meet difficulties.

What I understand:

  1. You have a non stationnary random process Y indexed by X such that the conditional expectation of Y given X=x[i] (if it is so, it's not the "mean" as you wrote it) is given by the formula CondExp(Y | X=x[i]) = a*x[i]+b, where a and b are known values.
     
  2. You assume this random process is a Gaussian randon process, then a realization of Y given X=x[i] (or [Y(omega) | X=x[i]] for short) is CondExp(Y | X=x[i]) + R(omega) where R is a gaussian random variable of 0 mean and standard deviation sigma.
     
  3. You want to generate a realization "GP" of this random process


If I'm right, and if X is some row vector of size N, just do

S := Statistics:-Sample(Normal(0, sigma), N);
GP := a *~ X +~ b +~ S

Maube you could try this:


A := [ some list on numerical values ]:
B := sort(A, output=permutation); # returns the position of the elements of A sorted in increasing values

For more info pease look to the 'sort' help page.
You'll see its possible to sort in descending order too

PS : once you have found, from some algorithm, the mosition of the highest element of a list A, the position of the second highest element of A is just the position of the highest element of the list A from which youy have removed the highest element.
 

To complete Preben's answer

Knowing that int( 1/(s*sqrt(2*Pi)) * exp(-1/2*(x/s)^2, x=-infinity..+infinity ) = 1, it's immediate to obtain, "by hand",  the value of this integral.

If You want to integrate over a bow in R3, the result is also obvious.
Just remark that  exp(-u^2/a) = exp(-(1/2) * (u /s))^2)  with s = sqrt(a/2).
Then  int( exp(-u^2/a), u=P..Q)  = s*sqrt(2*Pi) * int( 1/(s*sqrt(2*Pi)) * exp(-(1/2) * (u /s))^2) , u=P..Q) 
                                                   = s*sqrt(2*Pi) * (Probability(U <= Q) -Probability(U <= P) 
where U is a Gaussian random variable of expectation 0 and standard deviation s.

The result involves the error function 

I supppose the data you import are containded in X?
In your case X is column vector (a matrix with its number of column equal to 1).

Doing regression implies at least 2 quantities: the independant variable(s) and the dependant one.
However, in your case, it seems there is only one.
It may be that you have stacked the independant and dependant variables in the structure X?
In this case reshape X as desired, for instance by writting XY := Matrix(2, 1440, convert(X, list))^+ to obtain a 1440x2 matrix with the first column containing the independant variable (for instance) and the second column the dependant one.

Other points:

  • you use x in the expression of sigma and X elsewhere.
  • the quotes are unnecessary: just write gor instance mu := add(X[i], i=1..numelems(X)) / numelems(X)
  • in the computation of sigma, reuse mu: sigma := add((X[i]-mu)^2, i=1..numelems(X)) / numelems(X) *****
  • ***** last: the denominator in the epression of sigma must be numelems(X)-1 and not  numelems(X), unless 2880 represents the entire population
     

To conclude, unlesse you really want to use the explicit expressions for the mean and the standard deviation, it's better to use the correspondind Maple's function (just refer to the package Statistics (for instance))

with(Statistics);
mu:= Mean(X);    # or m := Statistics:-Mean(X) if you do not want to load the whole Statistics package
sigma := StandardDeviation(X)

 

 

Last point, when you Import yout data with ExcelTools:-Import (see the help page) you can specity the value an empty cell must contain.
Suppose this value is some number A which cannot be a "correct"number of your sample (for instance A := I, the imaginary root).
If X is the name of the matrix you inported, you can then extract the submatrix whose the components are real (or different from A in a more general way).

 

Could you please clarify your request?

By the way, I guees the integration range in the definition of M could be x=a..Q instead of x=0..Q ?

Could it be you are looking for that?

restart :
with(Statistics):

z := proc(S, c0, cs, F)
local a, b, M, N, ECost, DCost, f:
f := unapply(F, x):
a:=min(S):
b:=max(s)
M := int((Q-x)*f(x), x=a..Q):
N := int((Q-x)*f(x), x=Q..b):
ECost:=c0*M+cs*N:
DCost := diff(ECost, Q):
return fsolve(DCost=0, Q)
end proc:

Sam := Sample(Uniform(10, 10), 10):

# usage

z(Sam, 20, 25, 1);   # the last argument means “the function f(x) is a constant”
z(Sam, 20, 25, Pi);   # same result whatever the constant is
z(Sam, 20, 25, x);   # f(x) = x
z(Sam, 20, 25, exp(x));   # f(x) = exp(x)


# Naive coding

R := 1000:

AllQs := Vector[row](R):
for r from 1 to R do
   Sam := Sample(Uniform(10, 10), 10):
   AllQs[r] := z(Sam, 20, 25, 1);
end do:

# more elaborate coding

Sam := Sample(Uniform(10, 10), [R, 10]):   # a  matrix of R samples

AllQs := map(r -> z( entries(Sam[r, ..], nolist), 20, 25, 1);
 

Have you try this

restart:
roll := rand(0.0 .. 1.0):
N := 100:  # for example
s := [ seq(roll(), k=1..N) ]:
RandPerm := sort(s, output=permutation);


BTW : a random permutation of {1, 2, ...N} will ALWAYS return {1, 2, ..., N}. You must work with "ordered" objects such as lists, vectors, ..., not with sets

Add these two lines to your code to obtain the correct estimators.
(for more info look to mmcdara last post)

Unweighted_X_Sample := [ seq( X[n]$(round(Y[n]), n=1..numelems(X) ) ]:
DataSummary(X);

PS: works here because the weights are integers represented as floats.
In more general case (real weights) the procedure is easily generalizable (a very classical situation encountered by people who use MCMC methods for instance)

Hi,

 

restart:

X := [1, 2, 3, 4];

babyGrid := proc()
[ Grid:-MyNode(), X]
end proc;

Grid:-Launch(babyGrid, imports=['X'])


More info on the help pages

@LichengZhang

 

Sorry, I did not realize you wanted the TOTAL number of cycles.

As far as I know there is no function in Maple to do this, but you can obtain all the cycles by just "merging" the cycles of the cycle basis.
In a few words, CycleBasis returns a list of cycles, each of them being a list of vertices.
From this list you construct the set of edges that form the graph cycle.
Once done,  "merging in an adhoc way" some of  these sets will produce a new set that corresponds to a new cycle.

I think the code below returns the good result but I have no time to vericy this on a lot of graphs.

(if it's not the case I will correct it from my login  @mmcdara)

 

restart:

with(GraphTheory):
with(SpecialGraphs):

G := OctahedronGraph();
DrawGraph(G);
Cycles := CycleBasis(G);
N := numelems(Cycles);

GRAPHLN(undirected, unweighted, [1, 2, 3, 4, 5, 6], Array(%id = 18446744074375667950), `GRAPHLN/table/29`, 0)

 

 

[[1, 3, 2, 4], [1, 3, 2, 5], [1, 3, 2, 6], [1, 3, 5], [1, 3, 6], [1, 4, 5], [1, 4, 6]]

 

7

(1)

# close the list of vertices by adding the first vertex to its right

for n from 1 to N do
   C||n := [Cycles[n][], Cycles[n][1]];
end do;

[1, 3, 2, 4, 1]

 

[1, 3, 2, 5, 1]

 

[1, 3, 2, 6, 1]

 

[1, 3, 5, 1]

 

[1, 3, 6, 1]

 

[1, 4, 5, 1]

 

[1, 4, 6, 1]

(2)

# Create the sets of edges for each cycle

for n from 1 to N do
   S||n := map(m -> {(C||n)[m], (C||n)[m+1]}, {$1..numelems(C||n)-1});
end do;

{{1, 3}, {1, 4}, {2, 3}, {2, 4}}

 

{{1, 3}, {1, 5}, {2, 3}, {2, 5}}

 

{{1, 3}, {1, 6}, {2, 3}, {2, 6}}

 

{{1, 3}, {1, 5}, {3, 5}}

 

{{1, 3}, {1, 6}, {3, 6}}

 

{{1, 4}, {1, 5}, {4, 5}}

 

{{1, 4}, {1, 6}, {4, 6}}

(3)

S0 := {seq(S||n, n=1..N)}

{{{1, 3}, {1, 5}, {3, 5}}, {{1, 3}, {1, 6}, {3, 6}}, {{1, 4}, {1, 5}, {4, 5}}, {{1, 4}, {1, 6}, {4, 6}}, {{1, 3}, {1, 4}, {2, 3}, {2, 4}}, {{1, 3}, {1, 5}, {2, 3}, {2, 5}}, {{1, 3}, {1, 6}, {2, 3}, {2, 6}}}

(4)

with(combinat, powerset):

PS := powerset(S0):
NS := numelems(PS);

128

(5)

# Do some surgery while avoiding the empty set and the basis sets

AllCycles := NULL:
for ns from N+2 to NS do
   AllCycles := AllCycles, `union`(seq(PS[ns][k], k=1..nops(PS[ns]))) minus  `intersect`(seq(PS[ns][k], k=1..nops(PS[ns])))
end do:
AllCycles := {AllCycles}:
NAC := numelems(AllCycles);

120

(6)

{seq(`union`(op~(AllCycles[n])), n=1..NAC)} union {convert~(Cycles, set)[]};
TotalNumberOfCycles := numelems(%);

{{1, 3, 5}, {1, 3, 6}, {1, 4, 5}, {1, 4, 6}, {2, 3, 5}, {2, 3, 6}, {1, 2, 3, 4}, {1, 2, 3, 5}, {1, 2, 3, 6}, {1, 2, 4, 5}, {1, 2, 4, 6}, {1, 2, 5, 6}, {1, 3, 4, 5}, {1, 3, 4, 6}, {1, 3, 5, 6}, {1, 4, 5, 6}, {1, 2, 3, 4, 5}, {1, 2, 3, 4, 6}, {1, 2, 3, 5, 6}, {1, 2, 4, 5, 6}, {1, 3, 4, 5, 6}, {1, 2, 3, 4, 5, 6}}

 

22

(7)

 


 

Download TotalNumberOfCycles.mw

 

The "median" is either the midpoint of a segment or either a statistical quantity. As you seem to be interested in computing the "median" for an arbitrary number of points, I guess that you refer to the statistical quantity.

In this, if you do not already know this, you will find in the "Statistics" package two procedures named Median and Mean that you can assemble in a single procedure

medianmean := proc(x::list)
uses Statistics:
return [x, Median(x), Mean(x)]
end proc:


Nevertheless you must be very carefull with the "median".
In Probability/Statistics it is defined as the quantile of order 0.5, which means roughly that 50% of the observations (= of the numbers in your list) are smaller than the median and 50% larger.
For small samples, e.g [1, 2, 3, 8], the estimation of the median is a complex problem.
For instance every point in ]2, 3[ verify the requirement that 50% of the points are smaller than it: so what is the median? Your code returns 2.5, which is an acceptable value, but not more acceptable than 2.0001.
For very large samples this problems still exists but is often less crucial.
A robust estimation of the median can be obtained through resampling procedures (for instance Bootstrap, a procedure the Statistics package).





If the mean is not ambiguous, there is always a problem to define the median.

 

 

sol := dsolve(SYS2, numeric, parameters=[n]):

# set the value of n ... n=1 for instance
N:=1:
sol(parameters=[N]):


# "initialize" the computation
tmax := 600.:
sol(tmax):


# plot the solution
odeplot(sol, [t, x(t)], t=0..tmax, title=cat("n=", N))
 

 

 

I suggest you to replace
Interpolation:-Interpolate([seq(x, x = -2 .. 2, .1)], [seq(func(x, j), x = -2 .. 2, .1)], method = cubic)
by
CurveFitting:-Spline([seq(x, x = -2 .. 2, .1)], [seq(func(x, j), x = -2 .. 2, .1)], x, degree=1)  # I tested only with degree=1

I think you can also remove the option 'numeric' in the integration.

I keep having network to network transfer problem and I can't provide you an "improved" version of your code.
But the modifications are very light and I guess you will have no problem to use them by yourself.

Hope it will work.

PS: this doesn't really answer your original question about the "failing" you've met with Interpolation:-Interpolate

@maple2015 

 

Very optimistic of you to say it would take at least 10 minutes !
(it took more than 1.5 hour to me)

 

A hint (sorry I can't provide you the Maple code for a network to network transfert problem)

  1. Comment your last "Determinant(KFF)" line and deplace it by
    DetKFF := z -> Determinant(evalf(subs(P=z, KFF)));
  2. To have an idea of what DetKFF(z) looks like do
    plot('DetKFF(z)', z=0.01..5)
    Do not forget the quotes !
  3. It seems (zoom on the plot) that there is a first root around 0.15.
    Execute
    fsolve('DetKFF(z)', z=0.1)
    You should get the answer  0.16999...

 

It'is impossible to hope finding a root of Determinant(KFF).
Since the computation of it is extremely time consuming (I killed the command before it ended).

Other hints:

  1. do not try to start fsolve from z=0 (singularity for z=0)
  2. if you have a close look to the plot near z=0 you will find numerical oscillations: they can avoided by increasing the value of "Digits" (to the detriment to the efficiency of plot and fsolve).



BTW: I did a little mistake in the numerical values

node 1 = (0,0)

node 2 = (0,1)

node 3 = (5,6)

node 4 = (6,6)

1 2 Page 1 of 2