MaplePrimes Questions

Hi, 

I want to simulate a decision rule (let's say a production control strategy to fix the ideas) based on the outcomes of a binomial random variable.

One stage of this simulation involves operations which correspond to the NaiveProc described below. This procedure is a straightforward and rather naive translation in Maple of these operations.
It is given here to help you to understand what is to be done.

This coding being extremely slow, I implemented the required operations differently.
Given the fact that I'm only concerned by binomial samplings with outcomes 0 and 1, I basically replaced these samplings by three operations  U:=... , Got01:=... , Got0:= ... .
Each new implementation is about 200 times faster than the naive one.

My question is: Can we do still better?
For information the same simulation realized with R runs in less than 15s when REP = 10^8 (more than 100 times faster than what I'm capable to do with Maple)

Thanks for your help.


 

restart:

with(Statistics):

# This is the naive coding of the problem
#
# It has the advantage to explain clearly (at least I hope so) what is to be done.

NaiveProc := proc(REP, a, b, N)
  local roll, rep, p, K, N0, N1:
  roll := rand(a .. b);
  N0   := 0:
  N1   := 0:
  for rep from 1 to REP do
    p := roll();
    K := Sample(Binomial(N, p), 1)[1]:
    if K=0 then
      N0 := N0+1:
    elif K=1 then
      N1 := N1+1:
    end if:
  end do:
  N0, N1:
end proc:
    
CodeTools:-Usage(NaiveProc(10^3, 0, 5e-4, 25000));

memory used=78.92MiB, alloc change=68.00MiB, cpu time=1.21s, real time=6.30s, gc time=24.80ms

 

93, 80

(1)

# As the previous code is very slow I implemented several variants based on the fact
# that I care only on values of K equal to 0 and 1 and that it is then useless to sample
# a Binomial random variable to do this.
#
# Each of them is much faster than the naive coding... but can we even faster?
# For comparison, the same simulation realized in R takes less than 15s to run with REP=10^8

pi := (n, p, k) -> binomial(n, k) * (1-p)^(n-k) * p^k;


f1 := proc(REP, a, b, N)
  local Q, P0, P1, P01, U, Got01, Got0, N01, N0, N1:
  Q     := Sample(Uniform(a, b), REP)^+:
  P0    := pi~(N, Q, 0):
  P1    := pi~(N, Q, 1):
  P01   := P0 + P1:
  U     := Sample(Uniform(0, 1), REP)^+:
  Got01 := Select(t -> is(t <= 0), U-P01):
  Got0  := Select(t -> is(t <= 0), U-P0 );
  N01   := numelems(Got01):
  N0    := numelems(Got0):
  N1    := N01-N0:
  N0, N1;
end proc:

CodeTools:-Usage(f1(10^5, 0, 5e-4, 25000));

print():

f2 := proc(REP, a, b, N)
  local Q, P0, P1, P01, U, Got01, Got0, N01, N0, N1:
  Q     := Sample(Uniform(a, b), REP)^+:
  P0    := pi~(N, Q, 0):
  P1    := pi~(N, Q, 1):
  P01   := P0 + P1:
  U     := Sample(Uniform(0, 1), REP)^+:
  Got01 := select[flatten](`<=`, U-P01, 0):
  Got0  := select[flatten](`<=`, U-P0 , 0):
  N01   := numelems(Got01):
  N0    := numelems(Got0):
  N1    := N01-N0:
  N0, N1;
end proc:

CodeTools:-Usage(f2(10^5, 0, 5e-4, 25000));

print():

f3 := proc(REP, a, b, N)
  local Q, P0, P1, P01, U, Got01, Got0, N01, N0, N1:
  Q     := Sample(Uniform(a, b), REP)^+:
  P0    := pi~(N, Q, 0):
  P1    := pi~(N, Q, 1):
  P01   := P0 + P1:
  U     := Sample(Uniform(0, 1), REP)^+:
  Got01 := Vector(REP, i -> `if`(U[i] <= P01[i], 1, 0)):
  Got0  := Vector(REP, i -> `if`(U[i] <= P0 [i], 1, 0)):
  N01   := add(Got01):
  N0    := add(Got0):
  N1    := N01-N0:
  N0, N1;
end proc:

CodeTools:-Usage(f3(10^5, 0, 5e-4, 25000));

proc (n, p, k) options operator, arrow; binomial(n, k)*(1-p)^(n-k)*p^k end proc

 

memory used=376.83MiB, alloc change=296.77MiB, cpu time=4.12s, real time=4.03s, gc time=430.38ms

 

8118, 7979

 

 

memory used=173.26MiB, alloc change=0 bytes, cpu time=2.30s, real time=1.98s, gc time=461.22ms

 

7993, 7793

 

 

memory used=154.19MiB, alloc change=0 bytes, cpu time=2.28s, real time=2.01s, gc time=432.40ms

 

7989, 7869

(2)



 

 


 

Download MakeItFaster.mw

When Maple converts a matrix or a vector to Latex it uses latex array, which  is fine, but what it does is insert an exatra \\ at the end of the last row, just before the \end{array} 

This causes problem for some software which reads the latex file. It causes the software to take this as empty row and inserts a blank row at the end.

Here is a small example

A:=Matrix([[1,2],[3,4]]);
Physics:-Latex(A);

              \left[\begin{array}{cc}1 & 2 \\3 & 4 \\\end{array}\right]

Notice the extra \\  just before the \end{array} In Latex, this tells latex compiler to insert newline basically.

When opening the latex file using external software, which reads Latex and converts it to GUI, it show this

And now I have to go edit the latex by hand removing the extra \\ (inside the editor I can do global search and replace)

A better outout from Maple's Latex command would be the following

\left[\begin{array}{cc}1 & 2 \\3 & 4 \end{array}\right]

Even though both are valid Latex, the above is better.  I was wondering how hard it is to do this change.

All what the conversion needs to do is just check if the row in the matrix (or Vector) is the last one, and omit the extra \\ at the end of the last row. That is all. 

I could ofcourse write code to try to parse this each time I call Latex and to remove \\ myself in the program, but it will be much easier if Maple automatically did this.

 

ps. if the Latex conversion source code could be made public, others could help make improvement to it for free, benefiting Maple and everyone using it. Maplesoft would still approve the changes before checkin, and retain the copy right of the source code ofcourse. 

 

Maple 2020.1

 

Hello again,   I'm trying to write a procedure to multiply quaternions.   We know that they can be represented as a vector of length 4.   See Wikipedia.   My try  -

qaa=[8,4,6,2] and qab =[2,,4,6,8]
mult:=proc(a,b)

end proc;

maybe needs some looping or data from a table.  

Regards,

Matt

Hello All, is there a command/instruction that can tell me the phase margin, gain margin automatically after I have a bode plot, please?

Units are a neverending story, here's a variant of a problem raised before.

How can i get the minimum of 2 values with Units[Standard] loaded, when one of them could be zero (this loosing its units)?

See also 

  • https://www.mapleprimes.com/questions/229640-Units-Disappearing
  • https://www.mapleprimes.com/questions/230165-Error-in-UnitsTestDimensions-Invalid
  • https://www.mapleprimes.com/questions/230064-Compare-Similar-Units

with(Units[Standard])

a := 15*Unit('m')

15*Units:-Unit(m)

(1)

b := 13000*Unit('mm')

13000*Units:-Unit(mm)

(2)

min(a, b)

13*Units:-Unit(m)

(3)

c := 0*Unit('m')

0

(4)

min(a, c)

Error, (in Units:-Standard:-min) units with incompatible dimensions found: {1, length}

 

``


 

Download UnitsStandardCompare.mw


 

restart;

M__h := 0.352e-1;

0.352e-1

 

0.34e-1

 

0.8354e-1

 

0.96e-2

 

.123

 

0.7258e-1

 

0.214e-1

 

0.219e-1

 

.123

 

.7902

 

.11

 

0.136e-3

 

0.5e-1

 

0.8910e-1

 

0.45e-1

 

.7

 

.7214

 

1.354

 

0.235e-1

(1)

pdes := [diff(B(t, x), t) = M__h-beta__1*B(t, x)*G(t, x)/N__h+beta__2*B(t, x)*G(t, x)/N__h-mu__h*B(t, x)+sigma__h*E(t, x)*(diff(B(t, x), x, x)), diff(C(t, x), t) = beta__1*B(t, x)*G(t, x)/N__h-u[1]*C(t, x)/(1+C(t, x))-mu__h*C(t, x)*(diff(C(t, x), x, x)), diff(DD(t, x), t) = beta__2*DD(t, x)*G(t, x)/N__h-u[1]*DD(t, x)/(1+DD(t, x))-mu__h*DD(t, x)-delta__1*DD(t, x)*(diff(DD(t, x), x, x)), diff(E(t, x), t) = u[1]*C(t, x)/(1+C(t, x))+u[1]*DD(t, x)/(1+DD(t, x))-(mu__h+sigma__h)*E(t, x)*(diff(E(t, x), x, x)), diff(F(t, x), t) = M__b-beta__3*F(t, x)*C(t, x)/N__b+beta__4*F(t, x)*DD(t, x)/N__b-mu__b*F(t, x)*(diff(F(t, x), x, x)), diff(G(t, x), t) = beta__3*F(t, x)*C(t, x)/N__b+beta__4*F(t, x)*DD(t, x)/N__b-mu__b*G(t, x)*(diff(G(t, x), x, x))];

[diff(B(t, x), t) = 0.352e-1-0.891056911e-1*B(t, x)*G(t, x)-0.96e-2*B(t, x)+0.8910e-1*E(t, x)*(diff(diff(B(t, x), x), x)), diff(C(t, x), t) = .6791869919*B(t, x)*G(t, x)-0.45e-1*C(t, x)/(1+C(t, x))-0.96e-2*C(t, x)*(diff(diff(C(t, x), x), x)), diff(DD(t, x), t) = .5900813008*DD(t, x)*G(t, x)-0.45e-1*DD(t, x)/(1+DD(t, x))-0.96e-2*DD(t, x)-0.235e-1*DD(t, x)*(diff(diff(DD(t, x), x), x)), diff(E(t, x), t) = 0.45e-1*C(t, x)/(1+C(t, x))+0.45e-1*DD(t, x)/(1+DD(t, x))-0.9870e-1*E(t, x)*(diff(diff(E(t, x), x), x)), diff(F(t, x), t) = .7214-.1739837398*F(t, x)*C(t, x)+.1780487805*F(t, x)*DD(t, x)-1.354*F(t, x)*(diff(diff(F(t, x), x), x)), diff(G(t, x), t) = .1739837398*F(t, x)*C(t, x)+.1780487805*F(t, x)*DD(t, x)-1.354*G(t, x)*(diff(diff(G(t, x), x), x))]

(2)

bcs := [(D[2](B))(t, 0) = 0, (D[2](B))(t, 1) = 0, (D[2](C))(t, 0) = 0, (D[2](C))(t, 1) = 0, (D[2](DD))(t, 0) = 0, (D[2](DD))(t, 1) = 0, (D[2](E))(t, 0) = 0, (D[2](E))(t, 1) = 0, (D[2](F))(t, 0) = 0, (D[2](F))(t, 1) = 0, (D[2](G))(t, 0) = 0, (D[2](G))(t, 1) = 0, B(0, x) = 100, C(0, x) = 70, DD(0, x) = 50, E(0, x) = 70, F(0, x) = 100, G(0, x) = 70]

[(D[2](B))(t, 0) = 0, (D[2](B))(t, 1) = 0, (D[2](C))(t, 0) = 0, (D[2](C))(t, 1) = 0, (D[2](DD))(t, 0) = 0, (D[2](DD))(t, 1) = 0, (D[2](E))(t, 0) = 0, (D[2](E))(t, 1) = 0, (D[2](F))(t, 0) = 0, (D[2](F))(t, 1) = 0, (D[2](G))(t, 0) = 0, (D[2](G))(t, 1) = 0, B(0, x) = .100, C(0, x) = .70, DD(0, x) = .50, E(0, x) = .70, F(0, x) = .100, G(0, x) = .70]

(3)

sol := pdsolve(pdes, bcs, numeric);

module () local INFO; export plot, plot3d, animate, value, settings; option `Copyright (c) 2001 by Waterloo Maple Inc. All rights reserved.`; end module

(4)

sol:-plot3d([B(t, x), C(t, x)], t = 0 .. 20, x = 0 .. 20)

Error, (in pdsolve/numeric/plot3d) unable to compute solution for t>HFloat(0.25):
Newton iteration is not converging

 

``


 

Download spatial_1.mw

I want to vary t from -15 to -7 and from 7 to 15 how to write the Explore command?
example : Explore(Fig(t), t=-15..-7 and t=7..15); which does't work.  Thank you.

The center of mass for the density function p(r,θ,φ) = ln(r^2 + θ^2 + φ^2 +1) over the solid sphere radius 2.

I'm not sure how to do this on Maple.  

Many thanks!

Hi friends!

I have the following problem.

I'm trying to generate a matrix given a certain lenght n and a polynomial g(x) of degree r as follows:

If n=10, g(x)=x^4+3x^3+x^2+2x+1 and r=degree(g(x))=4, then the resulting matrix will be

 

 

 

 

 

i.e. a matrix with 10 columns, n-r=6 lines, and whose inputs are the coefficients of the polynomial, as follows:

Any help will be very appreciated.

Thank you!

I'm doing som XML readin, and need a bit help in using the HasChild function properly.

The attachment shows that using HasChild on a xmltree works fine, but I can't get i working when using it on a xmldocument. The latter is usually the case when you get when you read in from a file.

HasChild.mw

Hello there, 

Would you please tell me how to re-write the 'PMSM_v_eq' as 'PMSM_flux_eq_desired'? My simple attempt was using the 'solve()' command, but it failed. 

Here is the worksheet:


 

restart;

PMSM_v_eq := V__alphabeta(t) = R__s * i__alphabeta(t) + L__s*diff(i__alphabeta(t), t) + diff(lambda__alphabeta(t), t);

V__alphabeta(t) = R__s*i__alphabeta(t)+L__s*(diff(i__alphabeta(t), t))+diff(lambda__alphabeta(t), t)

(1)

 

PMSM_flux_eq := solve(PMSM_v_eq, lambda__alphabeta(t));

Error, (in solve) cannot solve expressions with diff(lambda__alphabeta(t), t) for lambda__alphabeta(t)

 

PMSM_flux_eq_desired := lambda__alphabeta(t) =  int(V__alphabeta(t) - R__s * i__alphabeta(t), t) - L__s*i__alphabeta(t);

lambda__alphabeta(t) = int(V__alphabeta(t)-R__s*i__alphabeta(t), t)-L__s*i__alphabeta(t)

(2)

 


Thank you!

Download PMSM_eq.mw

Hello

I have no choice but use Grid:-Map and Grid:-Seq in my calculations due to the size of them.  Here is a very small example that is puzzling me (Perhaps I did something really silly and did not realize). 

ansa:=CodeTools:-Usage(Grid:-Map(w->CondswithOnesolutionTest(w,eqns,vars,newvars,tlim),conds5s)):

with the following result:

ansa:=set([{alpha[1, 1] = 0, alpha[1, 2] = 0, alpha[1, 3] = 0, alpha[1, 4] = 0, alpha[1, 5] = 0, alpha[1, 6] = 0, alpha[1, 8] = 0, alpha[1, 9] = 0, alpha[2, 0] = 0, alpha[2, 1] = 0, alpha[2, 2] = 0, alpha[2, 4] = 0, alpha[2, 5] = 0, alpha[2, 7] = 0, alpha[2, 8] = 0, alpha[2, 9] = 0, alpha[3, 0] = 0, alpha[3, 1] = 0, alpha[3, 2] = 0, alpha[3, 3] = 0, alpha[3, 4] = 0, alpha[3, 6] = 0, alpha[3, 7] = 0, alpha[3, 8] = 0, alpha[3, 9] = 0}, {}, {}, {}, {}, {}], [{}, {}, {}, {}, {}, {alpha[1, 1] = 0, alpha[1, 2] = 0, alpha[1, 3] = 0, alpha[1, 4] = 0, alpha[1, 5] = 0, alpha[1, 6] = 0, alpha[1, 8] = 0, alpha[1, 9] = 0, alpha[2, 0] = 0, alpha[2, 1] = 0, alpha[2, 2] = 0, alpha[2, 4] = 0, alpha[2, 5] = 0, alpha[2, 7] = 0, alpha[2, 8] = 0, alpha[2, 9] = 0, alpha[3, 0] = 0, alpha[3, 1] = 0, alpha[3, 2] = 0, alpha[3, 3] = 0, alpha[3, 5] = 0, alpha[3, 6] = 0, alpha[3, 7] = 0, alpha[3, 8] = 0, alpha[3, 9] = 0}], [{alpha[1, 1] = 0, alpha[1, 2] = 0, alpha[1, 3] = 0, alpha[1, 4] = 0, alpha[1, 5] = 0, alpha[1, 6] = 0, alpha[1, 8] = 0, alpha[1, 9] = 0, alpha[2, 0] = 0, alpha[2, 1] = 0, alpha[2, 2] = 0, alpha[2, 4] = 0, alpha[2, 5] = 0, alpha[2, 7] = 0, alpha[2, 8] = 0, alpha[2, 9] = 0, alpha[3, 0] = 0, alpha[3, 1] = 0, alpha[3, 3] = 0, alpha[3, 4] = 0, alpha[3, 5] = 0, alpha[3, 6] = 0, alpha[3, 7] = 0, alpha[3, 8] = 0, alpha[3, 9] = 0}, {}, {}, {}, {}, {}], [{}, {}, {}, {}, {}, {alpha[1, 1] = 0, alpha[1, 2] = 0, alpha[1, 3] = 0, alpha[1, 4] = 0, alpha[1, 5] = 0, alpha[1, 6] = 0, alpha[1, 8] = 0, alpha[1, 9] = 0, alpha[2, 0] = 0, alpha[2, 1] = 0, alpha[2, 2] = 0, alpha[2, 4] = 0, alpha[2, 5] = 0, alpha[2, 7] = 0, alpha[2, 8] = 0, alpha[2, 9] = 0, alpha[3, 0] = 0, alpha[3, 2] = 0, alpha[3, 3] = 0, alpha[3, 4] = 0, alpha[3, 5] = 0, alpha[3, 6] = 0, alpha[3, 7] = 0, alpha[3, 8] = 0, alpha[3, 9] = 0}])

The same thing but now using only map

ansb:=CodeTools:-Usage(map(w->CondswithOnesolutionTest(w,eqns,vars,newvars,tlim),conds5s)):
ansb:={[{}, {}, {}, {}, {}, {alpha[1, 1] = 0, alpha[1, 2] = 0, alpha[1, 3] = 0, alpha[1, 4] = 0, alpha[1, 5] = 0, alpha[1, 6] = 0, alpha[1, 8] = 0, alpha[1, 9] = 0, alpha[2, 0] = 0, alpha[2, 1] = 0, alpha[2, 2] = 0, alpha[2, 4] = 0, alpha[2, 5] = 0, alpha[2, 7] = 0, alpha[2, 8] = 0, alpha[2, 9] = 0, alpha[3, 0] = 0, alpha[3, 1] = 0, alpha[3, 2] = 0, alpha[3, 3] = 0, alpha[3, 5] = 0, alpha[3, 6] = 0, alpha[3, 7] = 0, alpha[3, 8] = 0, alpha[3, 9] = 0}], [{}, {}, {}, {}, {}, {alpha[1, 1] = 0, alpha[1, 2] = 0, alpha[1, 3] = 0, alpha[1, 4] = 0, alpha[1, 5] = 0, alpha[1, 6] = 0, alpha[1, 8] = 0, alpha[1, 9] = 0, alpha[2, 0] = 0, alpha[2, 1] = 0, alpha[2, 2] = 0, alpha[2, 4] = 0, alpha[2, 5] = 0, alpha[2, 7] = 0, alpha[2, 8] = 0, alpha[2, 9] = 0, alpha[3, 0] = 0, alpha[3, 2] = 0, alpha[3, 3] = 0, alpha[3, 4] = 0, alpha[3, 5] = 0, alpha[3, 6] = 0, alpha[3, 7] = 0, alpha[3, 8] = 0, alpha[3, 9] = 0}], [{alpha[1, 1] = 0, alpha[1, 2] = 0, alpha[1, 3] = 0, alpha[1, 4] = 0, alpha[1, 5] = 0, alpha[1, 6] = 0, alpha[1, 8] = 0, alpha[1, 9] = 0, alpha[2, 0] = 0, alpha[2, 1] = 0, alpha[2, 2] = 0, alpha[2, 4] = 0, alpha[2, 5] = 0, alpha[2, 7] = 0, alpha[2, 8] = 0, alpha[2, 9] = 0, alpha[3, 0] = 0, alpha[3, 1] = 0, alpha[3, 2] = 0, alpha[3, 3] = 0, alpha[3, 4] = 0, alpha[3, 6] = 0, alpha[3, 7] = 0, alpha[3, 8] = 0, alpha[3, 9] = 0}, {}, {}, {}, {}, {}], [{alpha[1, 1] = 0, alpha[1, 2] = 0, alpha[1, 3] = 0, alpha[1, 4] = 0, alpha[1, 5] = 0, alpha[1, 6] = 0, alpha[1, 8] = 0, alpha[1, 9] = 0, alpha[2, 0] = 0, alpha[2, 1] = 0, alpha[2, 2] = 0, alpha[2, 4] = 0, alpha[2, 5] = 0, alpha[2, 7] = 0, alpha[2, 8] = 0, alpha[2, 9] = 0, alpha[3, 0] = 0, alpha[3, 1] = 0, alpha[3, 3] = 0, alpha[3, 4] = 0, alpha[3, 5] = 0, alpha[3, 6] = 0, alpha[3, 7] = 0, alpha[3, 8] = 0, alpha[3, 9] = 0}, {}, {}, {}, {}, {}]}

(This is what I expected as the result).

 

Why did Grid:-Map add set to the answer?  What am I missing?  

 

Many thanks

 

mwe.mw

I have a system of ODEs, and I am using the DifferentialAlgebra Package. I am getting this error message: Error, (in DifferentialAlgebra:-RosenfeldGroebner) unexpected occurrence of the non-rational constants {-I, I, 2*I} in the given input. What should I do?

Here is the problem. I start Maple 2020 on windows 10. Run a script which takes 1-2 days to complete. 

During this time, I can't use that Maple at all, since it is busy. 

I could start Maple 2019, and that runs as completely separate process. But I want to use Maple 2020 since some things in my scripts do not work on Maple 2019 that work on Maple 2020.

If I start a new instance of Maple 2020, by doing Start->Maple 2020. it does seem to start it OK, but I noticed it seems to be somehow still connected to the one running somehow.  May be they are sharing the same interface?

I can use the new instance now and open new worksheet and use it. But it seems to become very slow, as if it is sharing something with the other Maple 2020 running the long script which uses lots of resources. It is not RAM issue, I have 64 GB RAM, and there is plenty of free RAM left. 

When I close the new Maple 2020 workseet I started, I get a message asking if I want to save the worksheet that I have open from the earlier instance which is still running ! 

I say no ofcourse, as I do not want to terminate that instance, I want to keep it running until the script is completed.

My question is: Could someone may be explain exactly what happens when one starts new Maple 2020, while one is allready running? Why it seems they are sharing either the interface or something else.  How to start completely separate Maple 2020 instance on same PC while one is allready running?

With Mathematica, this issue does not happen. I can start two instances of same version on same PC, and there is nothing shared between them at all.  This does not seem to be the case with Maple.

Maple 2020.1 on windows 10.

 

sometimes solve returns solution of the form

restart;

#eqs:=.....
#sol:=solve(eqs,{v[1],v[2],v[3]});

sol:={v[1]=t,v[2]=3/2*t,v[3]=v[3]};

And wanted to remove all those that represnt arbitrary solution, which is v[3]=v[3] above.

I could do this using

remove(x->lhs(x)=rhs(x),sol);

which gives

{v[1]=t,v[2]=3/2*t};

But as an excersise, I could not figure how to do the same using subsindent (where I wanted to replace v[3]=v[3] with {} or NULL,. and also using applyrule.

Is it possible to do the same as above but using subsindent and applyrule (which is similar to patmatch)? 

First 484 485 486 487 488 489 490 Last Page 486 of 2427