tomleslie

9290 Reputation

17 Badges

11 years, 215 days

MaplePrimes Activity


These are answers submitted by tomleslie

well I don't read/speak Dutch, but the simplest way to "transform" plots is with the plottools:-transform() command. I guess the clue is in the name?

Maybe the attached provides what you want?

  restart;
  expr:=piecewise( x<=-2, undefined,
                   x<=0,  x^2+2*x,
                   x<=2,  -x^2+2*x,
                   undefined
                 );
  with(plottools):
  f1:=transform( (x,y)->[-x, y]):
  f2:=transform( (x,y)->[-2*x, y]):
  f3:=transform( (x,y)->[x-1, y]):
  f4:=transform( (x,y)->[2*x, y]):
  p1:=plot( expr, x=-2..2):
  plots:-display( [p1, f1(p1), f2(p1), f3(p1), f4(p1)],
                  color=[red, green, blue, cyan, black]
                )

expr := piecewise(x <= -2, undefined, x <= 0, x^2+2*x, x <= 2, -x^2+2*x, undefined)

 

 

 

Download transfPlt.mw

 

 

The DirectSearch package (available for free, from the Maple Application centre, thinks that there are 95 solutions none of which have a residual greater than ~10*(-32)

See the attached


 

restart;

[5*(0.2477538002e-5*a^2*b+0.1238769001e-5*b^3-0.1761804802e-9*b)*sin(Upsilon)+0.4081436800e-8*a-0.1717419154e-4*a^3, ((0.2642707202e-9*a*sigma+0.3320825870e-7*a+1.0*(-0.3716307003e-5*a*sigma-0.4669911380e-3*a)*b^2-0.1858153501e-5*a^3*sigma-0.2334955689e-3*a^3)*sin(Upsilon)+0.1280000000e-11*(5.026400000*a*sigma+631.6174240*a)*cos(Upsilon)+0.1280000000e-11*(1.352400000*sigma^2+339.8851680*sigma+21354.98511)*b+0.1280000000e-11*(-8.114400000*sigma^4-4078.62201638042*sigma^3-768779.463909262*sigma^2-64403218.29*sigma-2023227102.)*b^3)*csgn(sigma+125.6600000)-0.7680000000e-12*(sigma+125.6600000)^2*b, ((0.1238769001e-5*sigma+0.1556637127e-3)*b^4+((0.2477538002e-5*sigma-0.1556637126e-3)*a^2-0.1761804802e-9*sigma-0.2213883913e-7)*b^2-0.2334955690e-3*a^4+0.3320825870e-7*a^2)*cos(Upsilon)-0.8084703030e-9*a^2*sin(Upsilon)+(-0.7783185636e-4-0.2477538002e-5*sigma)*a*b^3+((0.3113274252e-3-0.1238769001e-5*sigma)*a^3+(-0.1106941957e-7+0.1761804802e-9*sigma)*a)*b+(0.1286758400e-8*sigma^2+0.1616940605e-6*sigma)*a*b]

proc (varr1) local vars, P; vars := [a, b, Upsilon]; eval(sys, sigma = varr1); fsolve(%); P := subs(%, vars) end proc

[0.4542414455e-3, -0.2720802421e-2, -27.31676601]

ans := DirectSearch:-SolveEquations(eval(sys, sigma = .5), AllSolutions = true, evaluationlimit = 100000)

_rtable[36893488148059668596]

(1)

 

380

(2)

 


 

Download DS.mw

 

  1. why you *want* to convert to a first order system (Maple will do this for you)
  2. precisley why you want to use a "shooting" method - Normally this would only be required if you wanted to adjust one of the parameters or initial conditions in order to achieve a specific outcome, for example maybe you want to find a value u0, in the initial condition u(0)=u0, in order that one obtains u(2)=3. Since you do not specify what you want to adjust, or the objective you wish to achieve, I can't help much with this

Otherwise the attached will solve your original ODE and pprint out a graph and table of results ofr any supple values of beta and u0


 

  restart;
  ode := diff(u(y), y, y) + beta*diff(diff(u(y), y)^2, y) = 2:
  bcs:= u(0)=u0, D(u)(0)=0:
  res:=dsolve([ode,bcs], numeric, parameters=[beta, u0]):

  getRes:=proc( p1, p2)
                local j:
                uses plots:
                res(parameters = [p1, p2]):
                printf
                ( "\n\n          beta=%12.8f u0=%12.8f\n",
                  p1,
                  p2
               ):     
                printf
                ( "%12a%18a%20a\n",
                  y,
                  u(y),
                  diff(u(y),y)
                ):
                seq
                ( printf
                  ( "%16.8f%16.8f%16.8f\n",
                    eval( y, res(j)),
                    eval(u(y),res(j)),
                    eval(diff(u(y),y), res(j))
                  ),
                  j=0..2, 0.1
                );
                odeplot
                ( res,
                  [ [y, u(y)],
                    [y, diff(u(y), y)]
                  ],
                  y=0..2,
                  color=[red, blue],
                  legend=[typeset(u(y)), typeset(diff(u(y),y))]
                );
  end proc:

  getRes(0.5, -5);
 



          beta=  0.50000000 u0= -5.00000000
           y              u(y)        diff(u(y),y)
      0.00000000     -5.00000000      0.00000000
      0.10000000     -4.99058294      0.18321599
      0.20000000     -4.96417443      0.34164083
      0.30000000     -4.92281209      0.48323974
      0.40000000     -4.86793763      0.61245159
      0.50000000     -4.80064122      0.73205089
      0.60000000     -4.72178493      0.84390900
      0.70000000     -4.63207269      0.94935888
      0.80000000     -4.53209347      1.04939024
      0.90000000     -4.42234981      1.14476112
      1.00000000     -4.30327656      1.23606805
      1.10000000     -4.17525560      1.32379005
      1.20000000     -4.03862490      1.40831907
      1.30000000     -3.89368730      1.48997992
      1.40000000     -3.74071539      1.56904670
      1.50000000     -3.57995657      1.64575142
      1.60000000     -3.41163713      1.72029411
      1.70000000     -3.23596412      1.79284821
      1.80000000     -3.05312861      1.86356440
      1.90000000     -2.86330794      1.93257567
      2.00000000     -2.66666650      2.00000005

 

 

 


 

Download odeProb.mw

 

for your second problem, I would prefer the graph shown in the attached: this provides the same information as that give by Acer's solution, it is just represented in a a different way.

  plot( [2^x, 1/2*x^2+2], x=-10..10, color=[red,blue]);
  fsolve( 2^x=1/2*x^2+2);

 

2.000000000

(1)

  plots:-complexplot( (-2)^x, x=-10..10, scaling=constrained);

 

 

Download dutch2.mw

 

perhaps? For some reason the final plot does not appear below, but I assure you that it does display in the worksheet!
 

restart

"x(t):= 2* exp(-t^(2));"

proc (t) options operator, arrow, function_assign; 2*exp(-t^2) end proc

(1)

p1 := plot(x(t), t = -3 .. 3, color = green, thickness = 4, size = [400, "square"])

p2 := plot(x(t), t = -3 .. 3, linestyle = "spacedash", color = blue, size = [400, "square"])

``

M := 3

K := [seq(k, k = -M .. M)]

xk := [seq(x(K[k]), k = 1 .. 2*M+1)]

NULLp3 := DynamicSystems[DiscretePlot](K, xk, style = stem, size = [400, "square"])

plots:-display([p1, p2, p3])

 

 


 

Download plts2.mw

  1. Doesn't seem to be an issue when the compile=false option is used (the default is compile=true)
  2. When compile=true seq() works OK, but add() crashes Maple

See the attached - don't uncomment the last command unless you happy to crash Maple

  kernelopts(version);
  with(Iterator):
#
# Set compile=false - no issues.
# NB, by default compile=true
#
  P := CartesianProduct([1,2,3,4],[1,2,3,4], compile=false);
  seq(p[], p=P);
  add(p[], p=P);

`Maple 2021.1, X86 64 WINDOWS, May 19 2021, Build ID 1539851`

 

CartesianProduct([1, 2, 3, 4], [1, 2, 3, 4], 'compile' = false)

 

Vector[row](2, {(1) = 1, (2) = 1}), Vector[row](2, {(1) = 2, (2) = 1}), Vector[row](2, {(1) = 3, (2) = 1}), Vector[row](2, {(1) = 4, (2) = 1}), Vector[row](2, {(1) = 4, (2) = 2}), Vector[row](2, {(1) = 3, (2) = 2}), Vector[row](2, {(1) = 2, (2) = 2}), Vector[row](2, {(1) = 1, (2) = 2}), Vector[row](2, {(1) = 1, (2) = 3}), Vector[row](2, {(1) = 2, (2) = 3}), Vector[row](2, {(1) = 3, (2) = 3}), Vector[row](2, {(1) = 4, (2) = 3}), Vector[row](2, {(1) = 4, (2) = 4}), Vector[row](2, {(1) = 3, (2) = 4}), Vector[row](2, {(1) = 2, (2) = 4}), Vector[row](2, {(1) = 1, (2) = 4})

 

Array(%id = 36893488148070420340)

(1)

#
# Use default setting so compile=true.
# seq() is OK but add() will cause
# Maple to crash, so is commented out
# below
#
  Q := CartesianProduct([1,2,3,4],[1,2,3,4], compile=false);
  seq(q[], q=Q);
#  add(q[], q=Q);

CartesianProduct([1, 2, 3, 4], [1, 2, 3, 4], 'compile' = false)

 

Array(%id = 36893488148070416244), Array(%id = 36893488148070408188), Array(%id = 36893488148070408308), Array(%id = 36893488148070408428), Array(%id = 36893488148070408548), Array(%id = 36893488148070408668), Array(%id = 36893488148070408788), Array(%id = 36893488148070408908), Array(%id = 36893488148070409028), Array(%id = 36893488148070409148), Array(%id = 36893488148070409268), Array(%id = 36893488148070409388), Array(%id = 36893488148070409508), Array(%id = 36893488148070409628), Array(%id = 36893488148070409748), Array(%id = 36893488148070409868)

(2)

 


 

Download atest.mw

when I set up a simple "toy" example, provided that the appropriate options are supplied to the contourplot() and impliicitplot() command. See the attached.

I can only suggest thta you use the big green up-arrow in the Mapleprimes toolbar to upload the worksheet in which your problem occurs

  restart:

#
# Use contourplot(). Note that by setting contours=[0]
# the command should return (all) contours on which
# sin(x*y)=0
#
  plots:-contourplot( sin(x*y),
                      x = -Pi..Pi,
                      y = -Pi..Pi,
                      contours=[0],
                      color=blue,
                      gridlines=false,
                      axes=boxed
                    );

 

#
# Use implicitplot() on the same expression. Note
# that since one is plotting an expression rather
# than an equation, Maple will *assume* that one
# wants to plot solution curves of the equation
# sin(x*y)=0 - which ought to be the same as the
# contourplot above
#
  plots:-implicitplot( sin(x*y),
                       x = -Pi..Pi,
                       y = -Pi..Pi,
                       color=blue,
                       gridlines=false,
                       axes=boxed
                     );

 

 

 

Download impCon2.mw

an expression such as 2.0+0.0*I is of type complex(extended, numeric).

On the other hand the result of simplify( 2.0+0..0*I), zero) is of type float.

See the attached, as well as the help for the command fnormal() below, my highlighting

fnormal preserves numeric type and sign information as much as possible.  Thus, for example, fnormal(1e-20*I, 10) = 0.*I not 0. .  This ensures that branching behavior is generally not affected by fnormal.  Use simplify(expr, zero) to remove 0 real or imaginary parts of complex floating point numbers.

  a:=2.0+0.0*I;
  whattype(a);
  b:=simplify(a, zero);
  whattype(b);

2.0+0.*I

 

complex(extended_numeric)

 

2.0

 

float

(1)

 

Download zeroComp.mw
 

to generate angle values in the range -720..-30, 30..720 and then seek to restrict these to "principal" values (which by definition will be in the range 0..360). For your case, only values in the range 30..360 wil be generated - so why generate any other values?

In any case, one of the three options in the attached may be what you want

  restart;
  randomize():
 

  L1:= [ seq( modp(j,360)*2*Pi/360,
              j=-720..-30, 15
            ),
         seq( modp(j,360)*2*Pi/360,
              j=30..720, 15
            )
       ]:
  r1:= rand(1..numelems(L1)):
  Da:= L1[r1()];

(5/3)*Pi

(1)

  L2:= [ seq( j,
              j=-720..-30, 15
            ),
         seq( j,
              j=30..720, 15
            )
       ]:
  r2:= rand(1..numelems(L2)):
  Da:= modp( L2[r2()], 360)*2*Pi/360;

(2/3)*Pi

(2)

  L3:= [ seq( j*2*Pi/360,
              j=30..360, 15
            )
       ]:
  r3:= rand(1..numelems(L3)):
  Da:= L3[r3()];

(19/12)*Pi

(3)

 

Download ang.mw

precisely what you want, but maybe one of the options in the attached?

  restart:

  eq1:=convert(E*D(II)(x)*D[1](psi)(x,t)+E*II(x)*(D[1]@@2)(psi)(x,t)-G*A(x)*psi(x,t)+G*A(x)*D[1](v)(x,t), diff)

E*(diff(II(x), x))*(diff(psi(x, t), x))+E*II(x)*(diff(diff(psi(x, t), x), x))-G*A(x)*psi(x, t)+G*A(x)*(diff(v(x, t), x))

(1)

  eq2:=convert(-(G*D[1](A)(x)*psi(x,t)+G*A(x)*D[1](psi)(x,t))+(G*D[1](A)(x)*v(x,t))+G*A(x)*(D[1]@@2)(v)(x,t)-m(x)/G*D[2](v)(x,t), diff);

-G*(diff(A(x), x))*psi(x, t)-G*A(x)*(diff(psi(x, t), x))+G*(diff(A(x), x))*v(x, t)+G*A(x)*(diff(diff(v(x, t), x), x))-m(x)*(diff(v(x, t), t))/G

(2)

  simplify(eval(eq1, psi(x,t)=rhs(isolate(eq2, psi))(x,t)));

(E*G^2*A(x)*(diff(A(x), x))^2*II(x)*(diff(diff(diff(diff(v(x, t), x), x), x), x))+2*E*(diff(A(x), x))*(-(diff(diff(A(x), x), x))*A(x)*II(x)+(diff(A(x), x))*((diff(A(x), x))*II(x)+(1/2)*(diff(II(x), x))*A(x)))*G^2*(diff(diff(diff(v(x, t), x), x), x))-E*G^2*A(x)*(diff(A(x), x))^2*II(x)*(diff(diff(diff(psi(x, t), x), x), x))-E*(diff(A(x), x))^2*II(x)*(diff(diff(diff(v(x, t), t), x), x))*m(x)+E*(diff(A(x), x))*II(x)*(G^2*A(x)*(diff(psi(x, t), x))-G^2*A(x)*(diff(diff(v(x, t), x), x))+m(x)*(diff(v(x, t), t)))*(diff(diff(diff(A(x), x), x), x))+(2*E*A(x)*II(x)*(diff(diff(A(x), x), x))^2-E*(diff(A(x), x))*((diff(A(x), x))*II(x)+(diff(II(x), x))*A(x))*(diff(diff(A(x), x), x))+(diff(A(x), x))^2*(E*(diff(II(x), x)+II(x))*(diff(A(x), x))-G*A(x)^2))*G^2*(diff(diff(v(x, t), x), x))-2*E*(diff(A(x), x))*(-(diff(diff(A(x), x), x))*A(x)*II(x)+(diff(A(x), x))*((diff(A(x), x))*II(x)+(1/2)*(diff(II(x), x))*A(x)))*G^2*(diff(diff(psi(x, t), x), x))-2*(-(diff(diff(A(x), x), x))*m(x)*II(x)+(diff(A(x), x))*((diff(m(x), x))*II(x)+(1/2)*(diff(II(x), x))*m(x)))*E*(diff(A(x), x))*(diff(diff(v(x, t), t), x))-2*E*II(x)*(G^2*A(x)*(diff(psi(x, t), x))+m(x)*(diff(v(x, t), t)))*(diff(diff(A(x), x), x))^2+((2*(diff(m(x), x))*II(x)+(diff(II(x), x))*m(x))*(diff(v(x, t), t))+(diff(psi(x, t), x))*G^2*((diff(A(x), x))*II(x)+(diff(II(x), x))*A(x)))*E*(diff(A(x), x))*(diff(diff(A(x), x), x))+(-E*(diff(v(x, t), t))*II(x)*(diff(diff(m(x), x), x))+(-E*(diff(m(x), x))*(diff(II(x), x))+G*A(x)*m(x))*(diff(v(x, t), t))+((-E*(diff(A(x), x))*(diff(II(x), x))+G*A(x)^2)*(diff(psi(x, t), x))+(diff(A(x), x))*((E*(diff(II(x), x))+G*A(x))*(diff(v(x, t), x))-G*A(x)*v(x, t)))*G^2)*(diff(A(x), x))^2)/(G^2*(diff(A(x), x))^3)

(3)

  eval(eq1, dsolve(eq2, psi(x,t)));

E*(diff(II(x), x))*((G^2*A(x)*(diff(diff(v(x, t), x), x))+G^2*(diff(A(x), x))*v(x, t)-m(x)*(diff(v(x, t), t)))/(G^2*A(x))-(Int((G^2*A(x)*(diff(diff(v(x, t), x), x))+G^2*(diff(A(x), x))*v(x, t)-m(x)*(diff(v(x, t), t)))/G^2, x)+_F1(t))*(diff(A(x), x))/A(x)^2)+E*II(x)*((G^2*(diff(A(x), x))*(diff(diff(v(x, t), x), x))+G^2*A(x)*(diff(diff(diff(v(x, t), x), x), x))+G^2*(diff(diff(A(x), x), x))*v(x, t)+G^2*(diff(A(x), x))*(diff(v(x, t), x))-(diff(m(x), x))*(diff(v(x, t), t))-m(x)*(diff(diff(v(x, t), t), x)))/(G^2*A(x))-2*(G^2*A(x)*(diff(diff(v(x, t), x), x))+G^2*(diff(A(x), x))*v(x, t)-m(x)*(diff(v(x, t), t)))*(diff(A(x), x))/(G^2*A(x)^2)+2*(Int((G^2*A(x)*(diff(diff(v(x, t), x), x))+G^2*(diff(A(x), x))*v(x, t)-m(x)*(diff(v(x, t), t)))/G^2, x)+_F1(t))*(diff(A(x), x))^2/A(x)^3-(Int((G^2*A(x)*(diff(diff(v(x, t), x), x))+G^2*(diff(A(x), x))*v(x, t)-m(x)*(diff(v(x, t), t)))/G^2, x)+_F1(t))*(diff(diff(A(x), x), x))/A(x)^2)-G*(Int((G^2*A(x)*(diff(diff(v(x, t), x), x))+G^2*(diff(A(x), x))*v(x, t)-m(x)*(diff(v(x, t), t)))/G^2, x)+_F1(t))+G*A(x)*(diff(v(x, t), x))

(4)

 

Download oddDiff.mw

is shown in the attached

  restart;
#
# May (or may not) want to randomize the seed
#
  randomize():
#
# The list of numbers
#
  L:= [ seq( `if`( j<>0, j , NULL),
             j=-720..720, 15
           )
      ]:
#
# Define a random number generator based on
# the list index
#
  r:= rand( 1..numelems(L) ):
#
# Generate a number
#
  L[r()];

-645

(1)

 

Download rand.mw

where I have fixed a couple of typos, tidied things up a little and commented out the random number generation in order to make the results reproducible (BTW your method of random number generation seems unnecessarily complicated).
 

NULL

NULL

restart

with(plots)

Da := 600

600

(1)

Ra := ((1/360)*Da*2)*Pi

(10/3)*Pi

(2)

xa := 9*cos(Ra)

-9/2

(3)

ya := 9*sin(Ra)

-(9/2)*3^(1/2)

(4)

arr := Ra-.2*Ra/abs(Ra)

10.27197551

(5)

``

xb := (1.5*((1/7)*abs(arr)+1))*cos(arr)

-2.450471313

(6)

yb := (1.5*((1/7)*abs(arr)+1))*sin(arr)

-2.773735738

(7)

xc := -.5*sin(Ra)*Ra/abs(Ra)

.2500000000*3^(1/2)

(8)

yc := .5*cos(Ra)*Ra/abs(Ra)

-.2500000000

(9)

Mulberry := ColorTools[Color]([169, 60, 147])

_m696303232

(10)

sketch1 :=  display([ textplot( [10, 0, 'x'],
                                align = {above, right}
                              ),
                      textplot( [0, 10, 'y'],
                                align = {above, right}
                              ),
                      plot( [t, 0, t = -10 .. 10],
                            thickness = 1,
                            color = black,
                            linestyle = solid
                          ),
                      plot( [0, t, t = -10 .. 10],
                            thickness = 1,
                            color = black,
                            linestyle = solid
                          ),
                      arrow( [xa, ya],
                             width = 0.1,
                             head_width = 0.5,
                             head_length = 0.5,
                             color = black
                           ),
                      plot( [ 1.5*(abs(t)/7 + 1)*cos(t),
                              1.5*(abs(t)/7 + 1)*sin(t),
                              t = 0 .. arr
                            ],
                            thickness = 3,
                            color = Mulberry,
                            linestyle = solid,
                            legend = beta,
                            legendstyle = [location = bottom]
                          ),
                      arrow( [xb, yb],
                             [xc, yc],
                             head_width = 0.5,
                             head_length = 0.5,
                             color = Mulberry
                           )
                    ],
                    axes = none,
                    labels = [``, ``],
                    view = [-10 .. 10, -10 .. 10]
                  );

 

NULL


 

Download aPlot2.mw

 

is wrong with using the command VariationalCalculus:-EulerLagrange(.....) ?

In your particular example

  1. I'd probablly use type/posint, rather than type/positive, because the latter doesn't care whether the argument is an integer or not
  2. the test you have written "works" but seems unnecessarliy complicated, and it seems you want it as an appllicable function, so I'd probably write it as in the attached

You can find out what Maple "knows" about an expression (including simple numbers) by using the op(0..) command as in the attached. I should warn you that AKAIK the existence of a zeroth operand is not "documented", but it generally "works".

See thet attached

  restart;
  f:= s -> type(s, posint)
           and
           type(log2(s), posint):
  f~([0, 4, -2]);

[false, true, false]

(1)

#
# An aside, what are the 'operands' associated
# with "simple" numbers in Maple
#
  op([0..], 2);
  op([0..], 2.0);
  op([0..], evalf(Pi));

Integer, 2

 

Float, 20, -1

 

Float, 3141592654, -9

(2)

 

Download posType.mw


although I think you might want to add some prettitying options to the "arrows" (PS the gridlines don't appear in the Maple worksheet, so the plot *looks* a bit cleaner

   restart;
   potfeld := (x, y) -> 1/sqrt(y^2+x^2):

   for i to 5 do
       for j to 3 do
           if   `and`( sqrt(i^2+j^2) <> 0,
                       sqrt((i-2)^2+j^2) <> 0,
                       evalf(potfeld(i, j)) < 3
                     )
           then P[i, j]:= plots:-arrow(`<,>`(i, j),
                                       `<,>`( (D[1](potfeld))(i, j),
                                              (D[2](potfeld))(i, j)
                                            )
                                      );
           else Pfeil[i, j]:= plots:-arrow( `<,>`(i, j),
                                            `<,>`(.1, .1)
                                          )
           end if;
       end do;
   end do;
   Pseq:= [ seq
            ( seq
              (  P[k, l],
                 k = 1 .. 5
              ),
              l = 1 .. 3
            )
          ]:
    plots:-display( Pseq,
                    view = [1 .. 5, 1 .. 3],
                    scaling = constrained
                  );

 

 

 


 

Download pfield.mw

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