tomleslie

11318 Reputation

18 Badges

12 years, 100 days

MaplePrimes Activity


These are answers submitted by tomleslie

as in the attached.
This particular command is very useful when dealing with commands which do not handle units

  restart:
  with(Units[Simple]):
  with(plots):

  a:=15*Unit('m'):
  b:=10*Unit('m'):
  displayPoints:=pointplot([a,b]);
  displayText := textplot( [ Units:-Split(a, output = coefficient),
                             Units:-Split(b, output = coefficient),
                             "text"
                           ]
                         );

 

 

 


 

Download Tplot2.mw

There are two attached woksheets with identical code - one of them was run in Maple 2021, the other in Maple 2022. It is fairly obvious that certain commands in the GroupTheory() package have changed between these two releases. Something as simple as what is returrned a the constructor command - sometimes it is a description and sometimes it is a list of generators.

A possible workaround for your issue in Maple 2021 is to "identify" the group in the database to which your DirectProduct() output is isomorphic, and apply the DrawSubGroupLattice() command to the former - this *seems* to work

  restart;
  interface(version);
  with(GroupTheory):
  G := DirectProduct(QuaternionGroup(), CyclicGroup(3));
  DrawSubgroupLattice(G);

`Standard Worksheet Interface, Maple 2022.0, Windows 7, March 8 2022 Build ID 1599809`

 

_m732551872

 

 

  id:=IdentifySmallGroup(G);
  G2:=SmallGroup(id);
  DrawSubgroupLattice(G2);

24, 11

 

_m755210368

 

 

 

Download GT2022.mw

 

  restart;
  interface(version);
  with(GroupTheory):
  G := DirectProduct(QuaternionGroup(), CyclicGroup(3));
  DrawSubgroupLattice(G);

`Standard Worksheet Interface, Maple 2021.2, Windows 7, November 23 2021 Build ID 1576349`

 

_m656624416

 

Error, (in .) `undefined[1]` does not evaluate to a module

 

  id:=IdentifySmallGroup(G);
  G2:=SmallGroup(id);
  DrawSubgroupLattice(G2);

24, 11

 

_m859921600

 

 

 

Download GT2021.mw

The attached using Maple 2022 seems to work just fine. Please specify the Maple version you are using and upload a worksheet (big greeen up-arrow in the Mapleprimes toolbar) illustrating your problem

  restart;
  with(GroupTheory):
  G := DirectProduct(QuaternionGroup(), CyclicGroup(3)):
  DrawSubgroupLattice(G)

 

  IsNormal(CyclicGroup(2), DirectProduct(CyclicGroup(2), CyclicGroup(2)));

true

(1)

 

Download GT.mw

I agree with pretty much everything Carl has said on this subject.

It is particularly importatnt to determine what happens at the "final" stage. As defined the system is "conservative": there is no dissipation, so the energy injected into the system by the iniitial impulse will "rattle around"  forever,

So knowing this, what would one expect to happen? If I apply some kind of "physical intuition", then the original impulse wouldprobably propagate forwards (and maybe backwards) through the system, untl this "original impulse" became "smoothed out" and pretty much everything in the system just kept bouncing around at the same (probably low-level), forever.

There is nothing particularly difficult about setting up and solving this system for (say) 50 stages see the attached. In the attached (I think) I have set up the system so that the right-hand end is a "brick wall" - in other words the right-hand end of the right-hand springs never move - for an N-stage system one just has to specify that x[N+1](t)=y[[N+1](t)=0, whch I have done in the attached.

This still does not expalin the fundamental question of the plots you are trying to produce. The "oscillator index" is an integer function, so I do not understand how the graphs you show have a "continuous" x-axis. What the hell is an "oscillator index" of (say) 15.333333? An "oscillator index" of 15, I understand: an "oscillator index" of 16, I understand: but an "!oscillator index" of (say) 15.333333 - I have no idea what this means!

In the attached, I have used 50 stages, and plotted a representain of the" kinetic energy" in each stage (mainly to show that this can be done). I have not attempted to address the potential energy stored in the springs at each stage. Of the three springs in each stage, two of them couple between stages, so the potential energy should probably be distributed between the two coupled stages.

Anyhow - for what it is worth, and with all the provisos given above, the attached just shows that it is perfectly possible to solve 100 coupled differential equations in Maple an plot some functions of the results. I'm obviously plotting the "wrong" function - but that is because I cannot work out what it is you want to plot :-(

  restart:
  with(plots):
  N:= 50: # number of stage in the system
  params:= [ m=0.022, k__g=1467.27, k__c=2.48e9, k__e=92.21, v=0.06]:
  odesys:= [ m*diff(x[1](t), t$2)+k__g*x[1](t)=-k__c*(x[1](t)-x[2](t))^3+k__e(y[1](t)-x[1](t)),
             m*diff(y[1](t), t$2)+k__g*y[1](t)=-k__c*(y[1](t)-y[2](t))^3+k__e(x[1](t)-y[1](t)),
             seq
             ( [ m*diff(x[j](t), t$2)+k__g*x[j](t)=k__c*((x[j-1](t)-x[j](t))^3-(x[j](t)-x[j+1](t))^3)+k__e(y[j](t)-x[j](t)),
                 m*diff(y[j](t), t$2)+k__g*y[j](t)=k__c*((y[j-1](t)-y[j](t))^3-(y[j](t)-y[j+1](t))^3)+k__e(x[j](t)-y[j](t))
               ][],
               j=2..N
             ),
             x[N+1](t)=0, # what happens at the right-hand end of the system??????
             y[N+1](t)=0
           ]:
  conds:= [        x[1](0)=0, D(x[1])(0)=v, y[1](0)=0, D(y[1])(0)=0,
            seq( [ x[j](0)=0, D(x[j])(0)=0, y[j](0)=0, D(y[j])(0)=0][],
                 j=2..N
               )
          ]:
  sol:= dsolve( eval( [odesys[], conds[]], params), numeric, maxfun=0):

#
# Plot *some sort* of representation of kinetic energy at each stage
#
  odeplot( sol,
           [ seq( [t, eval( (m/2)*(diff(x[j](t),t)^2+diff(y[j](t),t)^2), params)],
                  j=1..N
                )
           ],
           t=0..1.5
         );

 

 

 

 

Download bigSys.mw

What do you have set under

Tools->Options->Display->Font anti-aliasing

mine is set to Enabled. (If you change this I think you have to use "Apply Globally", then shut down Maple and restart for any changes to take effect

(sometimes using guesswork to figure out your intention) I can come up with the attached - which at least computes a solution for  F[3],F[4],F[5],........

However the equation to be solved at each stage contains the undefined function M(k+1). Your original code gives no indication of what this function might be, and I can't guess.

In future please use the big green up-arrow in the Mapleprimes toolbar to upload your worksheet. It save resdponders here from a lot of tedious, error-prone, retyping.

  restart;
  F[0]:=0: F[1]:=0: F[2]:=A:
  for k from -1 to 10 do
      delta[k]:=0:
  od:

  for k from 0 to 10 do
      F[k+3]:= solve
               ( F[k+3]+1/(k+1)/(k+2)/(k+3)*(1/2*add(cos(beta)*delta[k-m-1]*(m+1)*(m+2)*F[m+2], m=0..k)
                 -
                 1/2*sin(beta)*add(F[k-m]*(m+1)*(m+2)*F[m+2], m=0..k)+M(k+1)*F[k+1]),
                 F[k+3]
               );
od;
       

0

 

-(1/24)*M(2)*A

 

(1/60)*sin(beta)*A^2

 

(1/2880)*M(4)*M(2)*A

 

-(1/720)*sin(beta)*M(2)*A^2-(1/12600)*M(5)*sin(beta)*A^2

 

(11/20160)*sin(beta)^2*A^3-(1/967680)*M(6)*M(4)*M(2)*A

 

(1/90720)*sin(beta)*M(4)*M(2)*A^2+(1/48384)*sin(beta)*M(2)^2*A^2+(1/362880)*M(7)*sin(beta)*M(2)*A^2+(1/6350400)*M(7)*M(5)*sin(beta)*A^2

 

-(1/17280)*sin(beta)^2*A^3*M(2)-(11/4536000)*M(5)*sin(beta)^2*A^3-(11/14515200)*M(8)*sin(beta)^2*A^3+(1/696729600)*M(8)*M(6)*M(4)*M(2)*A

 

(5/266112)*sin(beta)^3*A^4-(29/958003200)*sin(beta)*M(6)*M(4)*M(2)*A^2-(7/22809600)*sin(beta)*M(4)*M(2)^2*A^2-(1/89812800)*M(9)*sin(beta)*M(4)*M(2)*A^2-(1/47900160)*M(9)*sin(beta)*M(2)^2*A^2-(1/359251200)*M(9)*M(7)*sin(beta)*M(2)*A^2-(1/6286896000)*M(9)*M(7)*M(5)*sin(beta)*A^2

 

(401/958003200)*M(4)*M(2)*A^3*sin(beta)^2+(563/319334400)*sin(beta)^2*M(2)^2*A^3+(37/479001600)*M(7)*sin(beta)^2*M(2)*A^3+(37/8382528000)*M(7)*M(5)*sin(beta)^2*A^3+(1/14784000)*M(2)*A^3*M(5)*sin(beta)^2+(1/22809600)*M(10)*sin(beta)^2*A^3*M(2)+(1/544320000)*M(10)*M(5)*sin(beta)^2*A^3+(1/1741824000)*M(10)*M(8)*sin(beta)^2*A^3-(1/919683072000)*M(10)*M(8)*M(6)*M(4)*M(2)*A

 

-(5023/2075673600)*sin(beta)^3*A^4*M(2)-(173/1945944000)*M(5)*sin(beta)^3*A^4-(23/1132185600)*M(8)*sin(beta)^3*A^4+(23/597793996800)*sin(beta)*M(8)*M(6)*M(4)*M(2)*A^2+(17/19926466560)*sin(beta)*M(2)^2*A^2*M(6)*M(4)+(1/948879360)*sin(beta)*M(4)^2*M(2)^2*A^2-(5/456648192)*M(11)*sin(beta)^3*A^4+(29/1643933491200)*M(11)*sin(beta)*M(6)*M(4)*M(2)*A^2+(7/39141273600)*M(11)*sin(beta)*M(4)*M(2)^2*A^2+(1/154118764800)*M(11)*M(9)*sin(beta)*M(4)*M(2)*A^2+(1/82196674560)*M(11)*M(9)*sin(beta)*M(2)^2*A^2+(1/616475059200)*M(11)*M(9)*M(7)*sin(beta)*M(2)*A^2+(1/10788313536000)*M(11)*M(9)*M(7)*M(5)*sin(beta)*A^2

(1)

 

Download loopSolve.mw

In the final 'for' loop, with i=0, Maple can find a symbolic solution for your ODE system. When i=1, Maple cannot find a symbolic solution so the dsolve() command returns an empty set. When you try to assign quantities such as

f1[i](x):=rhs(P[1]):

since P is the empty set, P[1] does not exist, and you get the error message

Error, invalid subscript selector

At this point I owuld normally suggest attempting a numerical solutionn for your ODE system, but since you have boundary conditions defined at four different values of the independent variable, this will not be possible.
##########################################

When attempting to run the code you supplied in more recent versions of Maple, many typesetting errors of the form

Error, (in Typeset:-Tdisplay[true]) improper op or subscript selector

Error, unexpected result from Typesetting

are produced. Changing the Typesetting level from 'extended' (default for versions later than Maple 18)  to 'standard' (the default for Maple 18) removes these typesetting errors, but introduces some odd problem with boundary conditions which I can't track down.

############################################

In the attached code I ahve "cleaned up" a lot of the notation without changing the 'logic' of the problem. This code will run the same way in Maple 18 with typesetting=standard and Maple 2022 with typesetting=extended, and produce the same "invalid subscript selector" in both.

  restart:
  with(plots):
  PDEtools[declare](f1(x),f2(x),f3(x), t1(x),t2(x),t3(x)):

f1(x)*`will now be displayed as`*f1

 

f2(x)*`will now be displayed as`*f2

 

f3(x)*`will now be displayed as`*f3

 

t1(x)*`will now be displayed as`*t1

 

t2(x)*`will now be displayed as`*t2

 

t3(x)*`will now be displayed as`*t3

(1)

  N :=1:
  F1:= add(p^jj*f1[jj](x), jj=0..N):
  F2:= add(p^jj*f2[jj](x), jj=0..N):
  F3:= add(p^jj*f3[jj](x), jj=0..N):
  T1:= add(p^jj*t1[jj](x), jj=0..N):
  T2:= add(p^jj*t2[jj](x), jj=0..N):
  T3:= add(p^jj*t3[jj](x), jj=0..N):

  gr:=5: pa:=5: sa:=4: br:=0.1: A1:=1: A2:=2: A3:=1:
  Eq11:= (1-p)*((diff(F1, x$2)+gr*T1)+pa)+p*((diff(F1, x$2)+gr*T1)+pa):
  Eq12:= (1-p)*(diff(T1, x$2))+p*((diff(T1, x$2)+br*(diff(F1, x))*(diff(F1, x)))):
  Eq21:= (1-p)*((diff(F2, x$2)+gr*A1*A2*T2)-sa*sa*F2+pa*A1)+p*((diff(F2, x$2)+gr*A1*A2*T2)-sa*sa*F2+pa*A1):
  Eq22:= (1-p)*(diff(T2, x$2))+p*((diff(T2, x$2)+A1*A3*br*((diff(F1, x))*(diff(F1, x))+sa*sa*F2*F2))):
  Eq31:= (1-p)*br*((diff(F3, x$2)+gr*T3)+pa)+p*((diff(F3, x$2)+gr*T3)+pa):
  Eq32:= (1-p)*br*(diff(T3, x$2))+p*((diff(T3, x$2)+br*(diff(F3, x))*(diff(F3, x)))):

  for i from 0 to N+1 do
      equ1[i] := coeff(Eq11, p, i) = 0:
      equ2[i] := coeff(Eq12, p, i) = 0:
      equ3[i] := coeff(Eq21, p, i) = 0:
      equ4[i] := coeff(Eq22, p, i) = 0:
      equ5[i] := coeff(Eq31, p, i) = 0:
      equ6[i] := coeff(Eq32, p, i) = 0:
  end do:

  con1[0]:= f1[0](-1) = 0, f1[0](0) = f2[0](0), D(f1[0])(0) = D(f2[0])(0),
            f2[0](1) = f3[0](1), D(f2[0])(1) = D(f3[0])(1),f3[0](2)=0:
  con2[0]:= t1[0](-1) = 0, t1[0](0) = t2[0](0), D(t1[0])(0) = D(t2[0])(0),
            t2[0](1) = t3[0](1), D(t2[0])(1) = D(t3[0])(1), t3[0](2)=1:
  for h to N do
      con1[h]:= f1[h](-1) = 0,  f1[h](0) = f2[h](0), D(f1[h])(0) = D(f2[h])(0),
                f2[h](1) = f3[h](1), D(f2[h])(1) = D(f3[h])(1), f3[h](2)=0:
      con2[h]:= t1[h](-1) = 0, t1[h](0) = t2[h](0), D(t1[h])(0) = D(t2[h])(0),
                t2[h](1) = t3[h](1), D(t2[h])(1) = D(t3[h])(1), t3[h](2)=0:
  end do:

  for i from 0 to N do
      P:= dsolve( {con1[i], con2[i], equ1[i], equ2[i], equ3[i], equ4[i], equ5[i], equ6[i]},
                  {f1[i](x), t1[i](x), f2[i](x), t2[i](x), f3[i](x), t3[i](x)}
                ):
      f1[i](x):=rhs(P[1]):
      f2[i](x):=rhs(P[2]):
      f3[i](x):=rhs(P[3]):
      t1[i](x):=rhs(P[4]):
      t2[i](x):=rhs(P[5]):
      t3[i](x):=rhs(P[6]):
  end do;

{f1[0](x) = -(5/18)*x^3-(10/3)*x^2-(5/72)*(447*exp(-4)-715*exp(4)+404)*x/(9*exp(-4)-25*exp(4))-(5/72)*(51*exp(-4)+404+385*exp(4))/(9*exp(-4)-25*exp(4)), f2[0](x) = -(5/144)*exp(4*x)*(237*exp(-4)+505)/(9*exp(-4)-25*exp(4))-(5/144)*exp(-4*x)*(303+395*exp(4))/(9*exp(-4)-25*exp(4))+(5/24)*x+25/48, f3[0](x) = -(5/18)*x^3-(10/3)*x^2+(5/72)*(316*exp(-4)*exp(4)+1605*exp(-4)-3785*exp(4))*x/(9*exp(-4)-25*exp(4))-(5/36)*(316*exp(-4)*exp(4)+597*exp(-4)-985*exp(4))/(9*exp(-4)-25*exp(4)), t1[0](x) = (1/3)*x+1/3, t2[0](x) = (1/3)*x+1/3, t3[0](x) = (1/3)*x+1/3}

 

-(5/18)*x^3-(10/3)*x^2-(5/72)*(447*exp(-4)-715*exp(4)+404)*x/(9*exp(-4)-25*exp(4))-(5/72)*(51*exp(-4)+404+385*exp(4))/(9*exp(-4)-25*exp(4))

 

-(5/144)*exp(4*x)*(237*exp(-4)+505)/(9*exp(-4)-25*exp(4))-(5/144)*exp(-4*x)*(303+395*exp(4))/(9*exp(-4)-25*exp(4))+(5/24)*x+25/48

 

-(5/18)*x^3-(10/3)*x^2+(5/72)*(316*exp(-4)*exp(4)+1605*exp(-4)-3785*exp(4))*x/(9*exp(-4)-25*exp(4))-(5/36)*(316*exp(-4)*exp(4)+597*exp(-4)-985*exp(4))/(9*exp(-4)-25*exp(4))

 

(1/3)*x+1/3

 

(1/3)*x+1/3

 

(1/3)*x+1/3

 

"P := "

 

Error, invalid subscript selector

 

 

Download badODEsys.mw

 

In your first execution group, for k=0 in the for loop, you are asking to solve the ode system

sys := {diff(c[0](x), x, x)/Pr, c[0](0) - 1, D(b[0])(5) - 1, diff(b[0](x), x, x, x), b[0](0), c[0](5), D(b[0])(0)}

In other words, two differential equations with dependent variables c[0](x) and b[0](x): system is second-order in c[0](x) and third-order in b[0](x). Sice you have two boudary conditions on c[0](x) and three boundary conditions on b[0](x), then there is sufficient, consistent, information for a solution, and Maple finds one.

################################################################################

In your second execution group, for k=0 in the for loop, you are asking to solve the ode system

sys:={1, diff(c[0](x), x, x)}

Now Maple does not like the fact that you have included the simple number '1' as part of an ODE system. What is it supposed to represent? Maple has no idea (and neither do I).

Maybe you intended a boundary condition such as c[0](0)=1 (in which case, you are still short of one boundary condition), or something else?

the source of your problem is in the plot commands where you have

eval(EQ2, params1 union {C = k/20})

If you check the indeterminates in the above expression with

indets(eval(EQ2, params1 union {C = k/20}));

this will return

{k, u, xi, z}

Now, k, u, z are determined within the plotting command either as range variables or the index of the seq() command - however 'xi' is completely undetermined, which is why you get the error message

Error, (in plots/implicitplot) invalid input: the following extra unknowns were found in the input expression: {xi}

Overall this seems like a complicated way to do phase portraits - if that is your ultimate objective. Have you considered eithe the phaseportrait()  or DEplot() commands? The attached shows how to produce some phase portrtaits and solution curves using the former - but I'm not sure if these satisfy your objective


 

  restart:
  with(DEtools):
  with(plots):

#
# Define system
#
  f := (u, z) -> z:
  g := (u, z) -> -u^3 + 3*z*u + c2*u^2/c1 + c2*z/c1 - c*u/c1:
  sys := diff(u(xi), xi) = f(u(xi), z(xi)),
         diff(z(xi), xi) = g(u(xi), z(xi)):
#
# Generate the phase portrait for a specified list of c, c1, c2
# with a couple of solution curves for different initial
# conditions
#
  phaseportrait
  ( eval
    ( {sys},
      [c=0.1, c1=-0.2, c2=-1]
    ),
    [u(xi), z(xi)],
    -1..1,
    [ [u(0) = 0, z(0) =  0.01],
      [u(0) = 0, z(0) = -0.01]
    ],
    color=blue,
    linecolor=[red, green]
  );
#
# Or solution curves for u(xi) only
#
  phaseportrait
  ( eval
    ( {sys},
      [c=0.1, c1=-0.2, c2=-1]
    ),
    [u(xi), z(xi)],
    -1..1,
    [ [u(0) = 0, z(0) =  0.01],
      [u(0) = 0, z(0) = -0.01]
    ],
    scene=[xi, u(xi)],
    color=blue,
    linecolor=[red, green]
  );
#
# Or solution curves for z(xi) only
#
  phaseportrait
  ( eval
    ( {sys},
      [c=0.1, c1=-0.2, c2=-1]
    ),
    [u(xi), z(xi)],
    -1..1,
    [ [u(0) = 0, z(0) =  0.01],
      [u(0) = 0, z(0) = -0.01]
    ],
    scene=[xi, z(xi)],
    color=blue,
    linecolor=[red, green]
  );

 

 

 

#
# The two solution curves in the above figures could also
# be produced with
#
  res1:= dsolve
         ( eval
           ( [sys, u(0)=u0,z(0)=z0 ],
             [u0=0, z0=0.01, c=0.1, c1=-0.2, c2=-1.0]
           ),
           numeric
         ):
  res2:= dsolve
         ( eval
           ( [sys, u(0)=u0,z(0)=z0],
             [u0=0, z0=-0.01, c=0.1, c1=-0.2, c2=-1.0]
           ),
           numeric
         ):
#
# z(xi) against u(xi)
#
  display
  ( [ odeplot
      ( res1,
        [u(xi), z(xi)],
        xi=-1..1,
        color=red
      ),
      odeplot
      ( res2,
        [u(xi), z(xi)],
        xi=-1..1,
        color=green
      )
    ]
  );
#
# u(xi) against xi
#
  display
  ( [ odeplot
      ( res1,
        [xi, u(xi)],
        xi=-1..1,
        color=red
      ),
      odeplot
      ( res2,
        [xi, u(xi)],
        xi=-1..1,
        color=green
      )
    ]
  );  
#
# z(xi) against xi
#
  display
  ( [ odeplot
      ( res1,
        [xi, z(xi)],
        xi=-1..1,
        color=red
      ),
      odeplot
      ( res2,
        [xi, z(xi)],
        xi=-1..1,
        color=green
      )
    ]
  );             

 

 

 

 

Download phPor.mw

 

 

the definition of ODE contains at least one (and very probably two) typos. The attached shows what the ODE should(?) be, although I guessed the "fix" for the second typo

It doesn't make a huge amount of difference because these is no way you are going to solve a single ODE which contains two unknown functions (along with their derivatives).

ODE := cos(g(t))^2*(diff(T(t), t, t))-3*rho*(diff(g(t), t))*cos(g(t))*sin(g(t))*(diff(T(t), t))+(Omega^2*cos(g(t))^z-8*rho^2*(diff(g(t), t))^2*sin(g(t))^2+2*sin(g(t))*cos(g(t))*(diff(g(t), t, t))*rho+2*(cos(g(t))^2)(diff(g(t), t))^2*rho+2*rho*(diff(g(t), t))^2*sin(g(t))^2)*T(t) = 0

cos(g(t))^2*(diff(diff(T(t), t), t))-3*rho*(diff(g(t), t))*cos(g(t))*sin(g(t))*(diff(T(t), t))+(Omega^2*cos(g(t))^z-8*rho^2*(diff(g(t), t))^2*sin(g(t))^2+2*sin(g(t))*cos(g(t))*(diff(diff(g(t), t), t))*rho+2*(cos(g(t)))(diff(g(t), t))^4*rho+2*rho*(diff(g(t), t))^2*sin(g(t))^2)*T(t) = 0

(1)

dsolve(ODE, T(t))

T(t) = DESol({-(-Omega^2*cos(g(t))^z+8*rho^2*(diff(g(t), t))^2*sin(g(t))^2-2*sin(g(t))*cos(g(t))*(diff(diff(g(t), t), t))*rho-2*(cos(g(t)))(diff(g(t), t))^4*rho-2*rho*(diff(g(t), t))^2*sin(g(t))^2)*_Y(t)/cos(g(t))^2-3*rho*(diff(g(t), t))*sin(g(t))*(diff(_Y(t), t))/cos(g(t))+diff(diff(_Y(t), t), t)}, {_Y(t)})

(2)

"#  #` Check unknowns in ODE - hmm contains a "product term" must be a typo:`  #` also contains a term in `(cos(g(t)))^(z) - is this for real"? or is 'z' another typo?"  #  indets(ODE);"

{Omega, rho, t, z, cos(g(t))^z, T(t), cos(g(t)), diff(T(t), t), diff(diff(T(t), t), t), diff(diff(g(t), t), t), diff(g(t), t), g(t), sin(g(t)), (cos(g(t)))(diff(g(t), t))}

(3)

ODE := cos(g(t))^2*(diff(T(t), t, t))-3*rho*(diff(g(t), t))*cos(g(t))*sin(g(t))*(diff(T(t), t))+(Omega^2*cos(g(t))^2-8*rho^2*(diff(g(t), t))^2*sin(g(t))^2+2*sin(g(t))*cos(g(t))*(diff(g(t), t, t))*rho+2*cos(g(t))^2*(diff(g(t), t))^2*rho+2*rho*(diff(g(t), t))^2*sin(g(t))^2)*T(t) = 0; indets(ODE)

cos(g(t))^2*(diff(diff(T(t), t), t))-3*rho*(diff(g(t), t))*cos(g(t))*sin(g(t))*(diff(T(t), t))+(Omega^2*cos(g(t))^2-8*rho^2*(diff(g(t), t))^2*sin(g(t))^2+2*sin(g(t))*cos(g(t))*(diff(diff(g(t), t), t))*rho+2*cos(g(t))^2*(diff(g(t), t))^2*rho+2*rho*(diff(g(t), t))^2*sin(g(t))^2)*T(t) = 0

 

{Omega, rho, t, T(t), cos(g(t)), diff(T(t), t), diff(diff(T(t), t), t), diff(diff(g(t), t), t), diff(g(t), t), g(t), sin(g(t))}

(4)

dsolve(ODE, T(t))

T(t) = DESol({diff(diff(_Y(t), t), t)-3*rho*(diff(g(t), t))*sin(g(t))*(diff(_Y(t), t))/cos(g(t))+(8*cos(g(t))^2*(diff(g(t), t))^2*rho^2+2*sin(g(t))*cos(g(t))*(diff(diff(g(t), t), t))*rho-8*rho^2*(diff(g(t), t))^2+Omega^2*cos(g(t))^2+2*rho*(diff(g(t), t))^2)*_Y(t)/cos(g(t))^2}, {_Y(t)})

(5)

 

Download badODE.mw

you have provided formulae for all the data you wish to compute, I'm having difficulty figuring out what your problem actually is??

Is it covered by the attached?

  restart;
  a := 1: b := 5:
  doDiff:= proc( ipData )
                 local j, opData:
                 if   numelems(ipData) > 1
                 then opData:=< seq
                                ( ipData[j]-ipData[j-1],
                                  j=2..numelems(ipData)
                                )
                              >:
                      return opData, thisproc(opData);
                 fi;
           end proc:
  xVals:= < $ a..b >:
  yVals:= < 1/~xVals >:
  yVals, doDiff(yVals);

Matrix(5, 1, {(1, 1) = 1, (2, 1) = 1/2, (3, 1) = 1/3, (4, 1) = 1/4, (5, 1) = 1/5}), Matrix(4, 1, {(1, 1) = -1/2, (2, 1) = -1/6, (3, 1) = -1/12, (4, 1) = -1/20}), Matrix(3, 1, {(1, 1) = 1/3, (2, 1) = 1/12, (3, 1) = 1/30}), Matrix(2, 1, {(1, 1) = -1/4, (2, 1) = -1/20}), Matrix(1, 1, {(1, 1) = 1/5})

(1)

 

Download doDiffs.mw

 

because I missed a term when typing in your equation.

When I include this term I see the same behaviour as you do, see the attached.

In the attached I also evaluated (the square root of) the discriminant, since this determines the nature of the possible solutions. This expression is actually

4*sqrt(-(a__02*b__11 - a__11*b__02)^2)

In other words, it is of the form

sqrt(-y^2)

Now if one applies the "simplification" simplify(expr, sqrt, symbolic), to this last expression, Maple will return I*y. Note that the help for simplify(expr, sqrt, symbolic) states that

The purpose of this feature is to allow simplification of expressions in contexts where the sign has no meaning, such as when x is an algebraic indeterminate.

So I think that Maple is correct - but I am somewhat surprised that it performs the last "simplification"


 

  restart;
  eqn:=4*a__02^2+4*k*a__02*a__11+k^2*a__11^2+4*k*b__02*b__11+k^2*b__11^2=0;
  solve(eqn, [k])[];

a__11^2*k^2+b__11^2*k^2+4*a__02*a__11*k+4*b__02*b__11*k+4*a__02^2 = 0

 

[k = 2*(-a__02*a__11-b__02*b__11+(-a__02^2*b__11^2+2*a__02*a__11*b__02*b__11+b__02^2*b__11^2)^(1/2))/(a__11^2+b__11^2)], [k = -2*(a__02*a__11+b__02*b__11+(-a__02^2*b__11^2+2*a__02*a__11*b__02*b__11+b__02^2*b__11^2)^(1/2))/(a__11^2+b__11^2)]

(1)

#
# Actually when typing 'eqn' above, I missed the term 4*b__02^2!
# When I include this term - the weird solve() behaviour
# appears
#
  eqn:=lhs(eqn)+4*b__02^2=0;
  solve(eqn, [k])[];

a__11^2*k^2+b__11^2*k^2+4*a__02*a__11*k+4*b__02*b__11*k+4*a__02^2+4*b__02^2 = 0

 

[k = 2*(I*a__02*b__11-I*a__11*b__02-a__02*a__11-b__02*b__11)/(a__11^2+b__11^2)], [k = -2*(I*a__02*b__11-I*a__11*b__02+a__02*a__11+b__02*b__11)/(a__11^2+b__11^2)]

(2)

#
# Idle curiosity - collect on 'k', just to show that eqn is
# a "simple" quadratic - makes no difference to the solution
#
  collect(eqn, k);
  solve(%, [k])[];

(a__11^2+b__11^2)*k^2+(4*a__02*a__11+4*b__02*b__11)*k+4*a__02^2+4*b__02^2 = 0

 

[k = 2*(I*a__02*b__11-I*a__11*b__02-a__02*a__11-b__02*b__11)/(a__11^2+b__11^2)], [k = -2*(I*a__02*b__11-I*a__11*b__02+a__02*a__11+b__02*b__11)/(a__11^2+b__11^2)]

(3)

#
# Take a look at (the square root of) the discriminant
#
  sqrt(discrim(lhs(eqn), k));
  simplify(%, sqrt, symbolic);

4*(-(a__02*b__11-a__11*b__02)^2)^(1/2)

 

(4*I)*b__11*a__02-(4*I)*a__11*b__02

(4)

 simplify(sqrt(-y^2), sqrt, symbolic);

I*y

(5)

 

Download solk2.mw

 

is quadratic in 'k' and therefore is trivial to solve (see attached). Whether the solutions are real or imaginary depends of all the parameters in the system.

As a general rule you should use the big green up-arrow in the Mapleprimes toolbar to upload your worksheet - it saves anyone here having to retype your stuff, which is time-consuming and error-prone


 

  restart;
  eqn:=4*a__02^2+4*k*a__02*a__11+k^2*a__11^2+4*k*b__02*b__11+k^2*b__11^2=0;
  solve(eqn, [k])[];

a__11^2*k^2+b__11^2*k^2+4*a__02*a__11*k+4*b__02*b__11*k+4*a__02^2 = 0

 

[k = 2*(-a__02*a__11-b__02*b__11+(-a__02^2*b__11^2+2*a__02*a__11*b__02*b__11+b__02^2*b__11^2)^(1/2))/(a__11^2+b__11^2)], [k = -2*(a__02*a__11+b__02*b__11+(-a__02^2*b__11^2+2*a__02*a__11*b__02*b__11+b__02^2*b__11^2)^(1/2))/(a__11^2+b__11^2)]

(1)

 

 

Download solk.mw

produces five (different) plots on my desktop.

It is probably(?) "safer" to use the plottools:-exportplot() command, rather than the more "generic" "Export()" command

  for i from 1 to 5 do
    plottools:-exportplot
               ( cat
                 ( "C:/Users/TomLeslie/Desktop/Fig",
                   i,
                   ".jpg"
                 ),
                 plot(x*i,x=0..10),
                 "JPEG"
               );
  od:

 

Download expPlot.mw

that syrup is clever enough to handle piecewise or conditional expressions within the circuit description itself.

However Syrup returns equations or ODEs which can be manipulated in the same way as any other Maple equation or ODE, so coming up with best case/worst case scenarios *ought* not to be too difficult - the optimisation package may be useful here.

I suggest you post a worksheet (use the big green up-arrow in the Mapleprimes toolbar) showing a "typical" example of what you are doing now, so that responders here can examine and make suggestions for improvements.

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