tomleslie

5595 Reputation

17 Badges

10 years, 7 days

MaplePrimes Activity


These are answers submitted by tomleslie

the attached is what you want

If not you are going to have to explain your issue more clearly
 

restart

with(plots)

with(CurveFitting)

Digits := 2

with(LinearAlgebra)

g := x -> piecewise(0 < x and x <= 10, x, 10 < x and x <= 20, -x+20);

g := proc (x) options operator, arrow; piecewise(0 < x and x <= 10, x, 10 < x and x <= 20, -x+20) end proc

(1)

p0 := plot(g(x), x = 0 .. 20)

 

plottools:-getdata(p0); M := %[-1]; ExportMatrix("E:/p0.xls", M)

["curve", [0. .. 20., 0. .. 9.99982415678391945], _rtable[18446744074378701510]]

 

_rtable[18446744074378701510]

(2)

v := .6

alpha := 15

mu := v/(2*alpha)

beta := v^2/(4*alpha)

NULL

#
# The following plots x-M*exp(-mu*x) for x=1.
# Fot other values of 'x' change the final option
# to x=2,2,-1, whatever
#
  plots:-pointplot( eval~(x-M*~exp(-mu*x), x=1));

 

``

Download plotMAt.mw

Achieving the attached convinced me that all calculations should be done "unitless" - and (only if strictly necessary) should units be applied at the end!

I'm not even sure that it is "complete" becuase the units of the expression returned by dsolve() ought to be 'm', but the expression *appears* to be "unitless"

Anyhow, for what it is worth

restart;
with(Units):
assume(t>0);

Automatically loading the Units[Simple] subpackage

 

 

T := 60*Unit('s');
mass := 80*Unit('kg');
b__1 := 13*Unit('kg')/Unit('s');
g := -9.81*Unit('m')/Unit('s')/Unit('s');
b := piecewise( 0 <= t*Unit('s') and t*Unit('s') < T,
                b__1,
                T<= t*Unit('s') and t*Unit('s') < infinity*Unit('s'),
                7*b__1
              );

Acceleration := diff(h(t)*Unit('m'), t*Unit('s'), t*Unit('s')) = g - b*diff(h(t)*Unit('m'), t*Unit('s'))/mass;

IC := h(0)*Unit('m') = 4000*Unit('m'),
      D(h)(0)*Unit(('m')/('s')) = 0.1*Unit('m'/'s');

solution := dsolve({IC, Acceleration}, h(t));

T := 60*Unit('s')

 

mass := 80*Unit('kg')

 

b__1 := 13*Unit('kg'/'s')

 

g := -9.81*Unit('m'/'s'^2)

 

b := piecewise(t < 60, 13, 60 <= t, 91)*Unit('kg'/'s')

 

Acceleration := diff(h(t), t, t) = -9.81-(1/80)*piecewise(t < 60, 13, 60 <= t, 91)*(diff(h(t), t))

 

IC := h(0) = 4000, (D(h))(0) = .1

 

h(t) = piecewise(t < 60, -(3924/65)*t-(62888/169)*exp(-(13/80)*t)+738888/169, 60 <= t, -(3924/455)*t+10118760/8281-(8984/169)*exp(-(91/80)*t+117/2)-(53904/169)*exp(-39/4)+(376704/8281)*exp(-(91/80)*t+273/4))

(1)

fsolve( rhs(solution), t);
plot( rhs(solution),
      t = 0 .. 150,
      y = 0 .. 4100,
      axis[2] = [gridlines = 30],
      axis[1] = [gridlines = 20],
      labels = ["t [s]", "h(t) [m]"],
      labelfont = ["Times", 15],
      color = "red", size = [1000, 1000]
    )

141.6838339

 

 

 


 

Download unitProb.mw

as in the attached. On this site, the animation appears to be runniing at 10 frames/sec (the default rate), but if you download the worksheet, you can use the worksheet animation toolbar to

  1. set the update rate to 1 frame/sec
  2. decide whether you wnat the animation to play a singel cycle or continuously
  3. manuallly step frame-by-frame

  restart;
  rnge:=-1..1:
  fcns:=[x, x^2, 1-x]:
  plots:-display( [ seq
                    ( plot
                      ( fcns[j], x=rnge),
                      j=1..3
                    )
                  ],
                  insequence=true
                );

 

 

 

 

 

Download slowanim.mw

One can use brute force, because it is relatively trivial to code, then check to see how bad timings get as the number of odd inputs and even addends is increased.

On my machine, the brute force approach in the attached produced the following approximate timings

#    Nodd    Neven      real time
#
#    100     100          15  msecs
#    100     1000        15  msecs
#    100     10000      15  msecs
#   
#    1000    100          66  msecs
#    1000    1000        77  msecs
#    1000    10000      1.0 secs
#  
#    10000   100         0.5 secs
#    10000   1000       1.5 secs
#    10000   10000      16  secs     

Out of idle curiosity, with Nodd=10000 and Neven =10000, there were 4534 "successful" odd numbers and 463 "failures".

A brief examination of the "failures" up to these limits indicates that they always have 3 as a prime factor - coincidence???

NB It would not surprise me if there were a much more efficient way to implement this test

  restart;
#
# Define the upper limit for the odd numbers
# to be tested and the upper limit for the
# even addends
#
  Nodd:=10000:
  Neven:=10000:
#
# Procedure which does the work
#
  doWork:= proc( p::posint, N::posint );
                 local f:=NumberTheory:-PrimeFactors(p),
                       j:
                 for j from 2 by 2 to N do
                     if   isprime~(f+~j)={true}
                     then #
                          # Return a list containing the
                          # input odd number, its prime
                          # factors, the lowest even number
                          # which produces a a new set of primes,
                          # and this new set of primes
                          #
                            return [p, f, j, f+~j];
                     fi;
                 od;
                 return #
                        # No even addend found: return a list
                        # containing the input odd number and
                        # a zero
                        #
                          [p, f, 0];
           end proc:

#
# Check how long it takes to return the complete list
# of answers for the specific values of NO and NE
#
# Results on my machine were. NB all timings approximate
#
#    Nodd    Neven      real time
#
#    100     100        15  msecs
#    100     1000       15  msecs
#    100     10000      15  msecs
#   
#    1000    100        66  msecs
#    1000    1000       77  msecs
#    1000    10000      1.0 secs
#  
#    10000   100        0.5 secs
#    10000   1000       1.5 secs
#    10000   10000      16  secs     
#
  ans:=CodeTools:-Usage([seq( doWork(j, Neven), j=3..Nodd,2)]):

memory used=1.73GiB, alloc change=40.00MiB, cpu time=14.98s, real time=15.04s, gc time=577.20ms

 

#
# Check answers for a couple of the cases used by OP
# to illustrate problem
#
# Utility to facilitate lookup of the answer in the
# main result list for any supplied odd number
#
  getVal:=p->(p-1)/2:
  ans[ getVal(119) ];
  ans[ getVal(105) ];

[119, {7, 17}, 6, {13, 23}]

 

[105, {3, 5, 7}, 0]

(1)

#
# Split the original list into successes and failures.
# Out of idle curiosity, check the number of successes
# and failures
#
  success,failure:=selectremove(i->numelems(i)=4, ans):
  numelems(success);
  numelems(failure);
#
# Output the first few successes and failures
#
  success[1..100];
  failure[1..100];
#
# Check whether the failure cases *always* have '3' as
# a prime factor. This should output any odd number
# which "failed" to generate a new set of prime numbers,
# and itself doesn't have '3' as a prime factor.
#
# It outputs nothing  interesting? coincidence??
#
  seq( `if`( j[2][1]=3, NULL, j[1]), j in failure);

4534

 

465

 

[[3, {3}, 2, {5}], [5, {5}, 2, {7}], [7, {7}, 4, {11}], [9, {3}, 2, {5}], [11, {11}, 2, {13}], [13, {13}, 4, {17}], [15, {3, 5}, 2, {5, 7}], [17, {17}, 2, {19}], [19, {19}, 4, {23}], [21, {3, 7}, 4, {7, 11}], [23, {23}, 6, {29}], [25, {5}, 2, {7}], [27, {3}, 2, {5}], [29, {29}, 2, {31}], [31, {31}, 6, {37}], [33, {3, 11}, 2, {5, 13}], [35, {5, 7}, 6, {11, 13}], [37, {37}, 4, {41}], [39, {3, 13}, 4, {7, 17}], [41, {41}, 2, {43}], [43, {43}, 4, {47}], [45, {3, 5}, 2, {5, 7}], [47, {47}, 6, {53}], [49, {7}, 4, {11}], [51, {3, 17}, 2, {5, 19}], [53, {53}, 6, {59}], [55, {5, 11}, 2, {7, 13}], [57, {3, 19}, 4, {7, 23}], [59, {59}, 2, {61}], [61, {61}, 6, {67}], [63, {3, 7}, 4, {7, 11}], [65, {5, 13}, 6, {11, 19}], [67, {67}, 4, {71}], [69, {3, 23}, 8, {11, 31}], [71, {71}, 2, {73}], [73, {73}, 6, {79}], [75, {3, 5}, 2, {5, 7}], [77, {7, 11}, 6, {13, 17}], [79, {79}, 4, {83}], [81, {3}, 2, {5}], [83, {83}, 6, {89}], [85, {5, 17}, 2, {7, 19}], [87, {3, 29}, 2, {5, 31}], [89, {89}, 8, {97}], [91, {7, 13}, 4, {11, 17}], [93, {3, 31}, 10, {13, 41}], [95, {5, 19}, 12, {17, 31}], [97, {97}, 4, {101}], [99, {3, 11}, 2, {5, 13}], [101, {101}, 2, {103}], [103, {103}, 4, {107}], [107, {107}, 2, {109}], [109, {109}, 4, {113}], [111, {3, 37}, 4, {7, 41}], [113, {113}, 14, {127}], [115, {5, 23}, 6, {11, 29}], [117, {3, 13}, 4, {7, 17}], [119, {7, 17}, 6, {13, 23}], [121, {11}, 2, {13}], [123, {3, 41}, 2, {5, 43}], [125, {5}, 2, {7}], [127, {127}, 4, {131}], [129, {3, 43}, 4, {7, 47}], [131, {131}, 6, {137}], [133, {7, 19}, 4, {11, 23}], [135, {3, 5}, 2, {5, 7}], [137, {137}, 2, {139}], [139, {139}, 10, {149}], [141, {3, 47}, 14, {17, 61}], [143, {11, 13}, 6, {17, 19}], [145, {5, 29}, 2, {7, 31}], [147, {3, 7}, 4, {7, 11}], [149, {149}, 2, {151}], [151, {151}, 6, {157}], [153, {3, 17}, 2, {5, 19}], [155, {5, 31}, 6, {11, 37}], [157, {157}, 6, {163}], [159, {3, 53}, 8, {11, 61}], [161, {7, 23}, 6, {13, 29}], [163, {163}, 4, {167}], [165, {3, 5, 11}, 2, {5, 7, 13}], [167, {167}, 6, {173}], [169, {13}, 4, {17}], [171, {3, 19}, 4, {7, 23}], [173, {173}, 6, {179}], [175, {5, 7}, 6, {11, 13}], [177, {3, 59}, 2, {5, 61}], [179, {179}, 2, {181}], [181, {181}, 10, {191}], [183, {3, 61}, 10, {13, 71}], [185, {5, 37}, 6, {11, 43}], [187, {11, 17}, 2, {13, 19}], [189, {3, 7}, 4, {7, 11}], [191, {191}, 2, {193}], [193, {193}, 4, {197}], [197, {197}, 2, {199}], [199, {199}, 12, {211}], [201, {3, 67}, 4, {7, 71}], [203, {7, 29}, 12, {19, 41}], [205, {5, 41}, 2, {7, 43}]]

 

[[105, {3, 5, 7}, 0], [195, {3, 5, 13}, 0], [231, {3, 7, 11}, 0], [285, {3, 5, 19}, 0], [315, {3, 5, 7}, 0], [357, {3, 7, 17}, 0], [429, {3, 11, 13}, 0], [465, {3, 5, 31}, 0], [483, {3, 7, 23}, 0], [525, {3, 5, 7}, 0], [555, {3, 5, 37}, 0], [585, {3, 5, 13}, 0], [609, {3, 7, 29}, 0], [627, {3, 11, 19}, 0], [645, {3, 5, 43}, 0], [663, {3, 13, 17}, 0], [693, {3, 7, 11}, 0], [735, {3, 5, 7}, 0], [855, {3, 5, 19}, 0], [861, {3, 7, 41}, 0], [897, {3, 13, 23}, 0], [915, {3, 5, 61}, 0], [945, {3, 5, 7}, 0], [969, {3, 17, 19}, 0], [975, {3, 5, 13}, 0], [987, {3, 7, 47}, 0], [1005, {3, 5, 67}, 0], [1023, {3, 11, 31}, 0], [1071, {3, 7, 17}, 0], [1095, {3, 5, 73}, 0], [1113, {3, 7, 53}, 0], [1131, {3, 13, 29}, 0], [1155, {3, 5, 7, 11}, 0], [1185, {3, 5, 79}, 0], [1221, {3, 11, 37}, 0], [1239, {3, 7, 59}, 0], [1287, {3, 11, 13}, 0], [1311, {3, 19, 23}, 0], [1365, {3, 5, 7, 13}, 0], [1395, {3, 5, 31}, 0], [1419, {3, 11, 43}, 0], [1425, {3, 5, 19}, 0], [1449, {3, 7, 23}, 0], [1455, {3, 5, 97}, 0], [1491, {3, 7, 71}, 0], [1545, {3, 5, 103}, 0], [1575, {3, 5, 7}, 0], [1581, {3, 17, 31}, 0], [1599, {3, 13, 41}, 0], [1617, {3, 7, 11}, 0], [1635, {3, 5, 109}, 0], [1653, {3, 19, 29}, 0], [1665, {3, 5, 37}, 0], [1743, {3, 7, 83}, 0], [1755, {3, 5, 13}, 0], [1785, {3, 5, 7, 17}, 0], [1827, {3, 7, 29}, 0], [1833, {3, 13, 47}, 0], [1869, {3, 7, 89}, 0], [1881, {3, 11, 19}, 0], [1887, {3, 17, 37}, 0], [1905, {3, 5, 127}, 0], [1935, {3, 5, 43}, 0], [1989, {3, 13, 17}, 0], [1995, {3, 5, 7, 19}, 0], [2013, {3, 11, 61}, 0], [2067, {3, 13, 53}, 0], [2079, {3, 7, 11}, 0], [2085, {3, 5, 139}, 0], [2121, {3, 7, 101}, 0], [2139, {3, 23, 31}, 0], [2145, {3, 5, 11, 13}, 0], [2193, {3, 17, 43}, 0], [2205, {3, 5, 7}, 0], [2211, {3, 11, 67}, 0], [2247, {3, 7, 107}, 0], [2265, {3, 5, 151}, 0], [2301, {3, 13, 59}, 0], [2325, {3, 5, 31}, 0], [2337, {3, 19, 41}, 0], [2355, {3, 5, 157}, 0], [2373, {3, 7, 113}, 0], [2409, {3, 11, 73}, 0], [2415, {3, 5, 7, 23}, 0], [2445, {3, 5, 163}, 0], [2499, {3, 7, 17}, 0], [2535, {3, 5, 13}, 0], [2541, {3, 7, 11}, 0], [2553, {3, 23, 37}, 0], [2565, {3, 5, 19}, 0], [2583, {3, 7, 41}, 0], [2607, {3, 11, 79}, 0], [2625, {3, 5, 7}, 0], [2679, {3, 19, 47}, 0], [2691, {3, 13, 23}, 0], [2697, {3, 29, 31}, 0], [2715, {3, 5, 181}, 0], [2745, {3, 5, 61}, 0], [2751, {3, 7, 131}, 0], [2769, {3, 13, 71}, 0]]

(2)

 

Download prmeProb.mw

(I'm on a Windows machine, so I can't debug/reproduce fully)

  1. Open relevant help page with ?Physics,Geodesics
  2. In the help page menu system, ensure that "View->Display Examples with 2D math input" is unchecked. In other words  all the executable math in the help page should be displayed as 1-D input
  3. In the help page menu system, use View->Open Page as Worksheet. which does exactly what it says. It will open a worksheet called [Physics, Geodesics], which you can execute in the standard way
  4. Execute this worksheet one "group" at a time - What happens?
  5. If this "works" you can try steps 1-4 gain, except at step (2) above, ensure that "View->Display Examples with 2D math input" is checked. - Again, What happens?

the conditions f(0)=0 and f'(0)=0 do not impose any restrictions on the values of the parameters n[1], n[2], n[3].

See the attached

   restart;
#
# Define the function
#
  f:=x-> 0.8680555556e-1*x*(3.262720*x^5*n[2]-.63360*x^5*n[1]^2+2.69184*x^5*n[1]-.480*x^3*n[1]^2+2.5920*x^3*n[2]+5.760*n[1]*x+1.920*n[2]*x^2+0.96e-1*x^4*n[1]^3+3.04320*x^4*n[2]-.5280*x^4*n[1]^2-0.32e-1*x^5*n[1]^4+.2976*x^5*n[1]^3-.1184*x^5*n[2]^2+.5376*x^5*n[3]-0.576e-1*x^4*n[1]+.576*x^4*n[3]-6.3360*x^3-6.96960*x^4-7.349760*x^5+11.520-.192*x^5*n[1]*n[3]-.7104*x^4*n[1]*n[2]+.2528*x^5*n[1]^2*n[2]-1.81824*x^5*n[1]*n[2]):
#
# Rewrite the above equation in a somewhat "neater"
# form by collecting powers of 'x', just to make
# it easier to read/interpret
#
  collect(expand(f(x)), x);
#
# Evaluate f(x) at x=0
#
  eval(f(x), x=0);
#
# Differentiate the above expression wrt 'x' and
# evaluate at x=0
#
  D(f)(0);

(.2832222222*n[2]-0.5500000000e-1*n[1]^2+.2336666667*n[1]-0.2777777778e-2*n[1]^4+0.2583333333e-1*n[1]^3-0.1027777778e-1*n[2]^2+0.4666666667e-1*n[3]-.6380000000-0.1666666667e-1*n[1]*n[3]+0.2194444445e-1*n[1]^2*n[2]-.1578333333*n[1]*n[2])*x^6+(0.8333333334e-2*n[1]^3+.2641666667*n[2]-0.4583333334e-1*n[1]^2-0.5000000000e-2*n[1]+0.5000000000e-1*n[3]-.6050000000-0.6166666667e-1*n[1]*n[2])*x^5+(-0.4166666667e-1*n[1]^2+.2250000000*n[2]-.5500000000)*x^4+.1666666667*x^3*n[2]+.5000000000*n[1]*x^2+1.000000000*x

 

0.

 

1.000000000

(1)

 


 

Download poly.mw

which

  1. plots all dependent variables
  2. produces a matrix of dependent variable values for T=0..20

restart; _local(gamma)

M__h := .50; beta__o := 0.34e-1; beta__1 := 0.25e-1; mu__r := 0.4e-3; sigma := .7902; alpha := .11; psi := 0.136e-3; xi := 0.5e-1; gamma := .7; M__c := .636; mu__b := 0.5e-2

ODEs := diff(B(T), T) = M__h-beta__1*psi*B(T)*G(T)-mu__r*B(T), diff(C(T), T) = beta__1*psi*B(T)*G(T)-sigma*psi*beta__1*E(T)*DD(T)-(alpha+xi+mu__r)*C(T), diff(DD(T), T) = alpha*C(T)-(gamma+mu__r)*DD(T), diff(E(T), T) = gamma*DD(T)+sigma*psi*beta__1*E(T)*G(T)-mu__r*E(T), diff(F(T), T) = M__c-psi*beta__o*F(T)*C(T)-mu__b*F(T), diff(G(T), T) = psi*beta__o*F(T)*C(T)-mu__b*G(T); bcs := B(0) = .50, C(0) = .30, DD(0) = .21, E(0) = .14, F(0) = .70, G(0) = .45

ans := dsolve([ODEs, bcs], numeric); for j in indets([ODEs], function(name)) do plots:-odeplot(ans, [T, j], T = 0 .. 30, title = convert(j, string), titlefont = [tims, bold, 20], thickness = 2, color = red) end do

 

 

 

 

 

 

interface(rtablesize = [22, 10]); M := Matrix([`~`[lhs](ans(0)), seq(`~`[rhs](ans(j)), j = 0 .. 20)])

Matrix(%id = 18446744074595705910)

(1)

 

``


 

Download odeRes.mw

just because I like alternatives

  restart;
#
# Random matrix for test purposes
#
  A:=LinearAlgebra:-RandomMatrix(10):
#
# Sort columns of A in ascending order
#
  B:=Matrix([seq(sort(A[..,j]), j=1..10)], scan=rows);
#
# Define "finder" function
#
  g:=(mat, col, val)-> ListTools:-SelectFirst
                       ( x-> evalb( x>val ),
                         mat[..,col],
                         output=indices
                       ):
#
# Usage example
#
# In the matrix B, output the row index of the
# first entry in column 2 which is >50
#
  g(B, 2, 50);
#
# Returns NULL if no entry exceeds supplied value
#
  g(B, 2, 100);

Matrix(10, 10, {(1, 1) = -63, (1, 2) = -35, (1, 3) = -82, (1, 4) = -68, (1, 5) = -95, (1, 6) = -61, (1, 7) = -80, (1, 8) = -81, (1, 9) = -98, (1, 10) = -31, (2, 1) = -38, (2, 2) = -14, (2, 3) = -70, (2, 4) = -67, (2, 5) = -25, (2, 6) = -48, (2, 7) = -50, (2, 8) = -50, (2, 9) = -93, (2, 10) = -4, (3, 1) = -26, (3, 2) = 12, (3, 3) = -32, (3, 4) = -62, (3, 5) = -20, (3, 6) = -44, (3, 7) = -2, (3, 8) = -38, (3, 9) = -77, (3, 10) = 8, (4, 1) = -23, (4, 2) = 19, (4, 3) = -13, (4, 4) = -59, (4, 5) = 9, (4, 6) = 9, (4, 7) = 10, (4, 8) = -22, (4, 9) = -76, (4, 10) = 27, (5, 1) = -1, (5, 2) = 21, (5, 3) = -1, (5, 4) = -33, (5, 5) = 14, (5, 6) = 20, (5, 7) = 12, (5, 8) = -18, (5, 9) = -74, (5, 10) = 29, (6, 1) = 10, (6, 2) = 45, (6, 3) = 29, (6, 4) = 12, (6, 5) = 16, (6, 6) = 24, (6, 7) = 25, (6, 8) = -16, (6, 9) = -72, (6, 10) = 44, (7, 1) = 22, (7, 2) = 60, (7, 3) = 41, (7, 4) = 18, (7, 5) = 22, (7, 6) = 65, (7, 7) = 31, (7, 8) = -9, (7, 9) = -32, (7, 10) = 67, (8, 1) = 30, (8, 2) = 80, (8, 3) = 52, (8, 4) = 42, (8, 5) = 51, (8, 6) = 76, (8, 7) = 43, (8, 8) = 33, (8, 9) = -2, (8, 10) = 69, (9, 1) = 63, (9, 2) = 88, (9, 3) = 70, (9, 4) = 72, (9, 5) = 60, (9, 6) = 77, (9, 7) = 50, (9, 8) = 45, (9, 9) = 27, (9, 10) = 92, (10, 1) = 91, (10, 2) = 90, (10, 3) = 91, (10, 4) = 82, (10, 5) = 99, (10, 6) = 86, (10, 7) = 94, (10, 8) = 87, (10, 9) = 57, (10, 10) = 99})

 

7

(1)

 


 

Download searchMat.mw

is shown in the attached

restart;
f:=k->sum((r+1)*(r+2)*(r+3)*(r+4)*F(r+4)*(k-r+1)*F(k-r+1), r = 0 .. k) =c:
g:=h->isolate( h, indets(h, function)[-1]):
(g @ f)(0);
(g @ f)(1);
(g @ f)(2);

F(4) = (1/24)*c/F(1)

 

F(5) = (1/120)*(c-48*F(4)*F(2))/F(1)

 

F(6) = (1/360)*(c-72*F(4)*F(3)-240*F(5)*F(2))/F(1)

(1)

 

Download altSum.mw

If I strip out all the redundant stuff you neither want nor need the I am left wiith the worksheet shown in the attached.

There are two reasons why no solution can be obtained from this worksheet, both of which you will have to fix

  1. You have six ODES , but seven independent variables, namely A(T), B(T), C(T), G(T), R(T), S(T), I(T). Rouben has already pointed out that one of these dependent variables ( ie I(T) ) has a name which will cause problems. You need to work out if you really do have seven dependent variables: if you do, then don;'t call one of them I(T) - any other name will be fine: eg P(T), Q(T), Z(T), whatever, just not I(T). However if you do have seven dependent variables, then either
    1. One of these is going to have to be defined in terms of the other six, or
    2. You are going to need another equation
  2. The variable mu__r is used in a couple of these equations and is never defined

see the attached

  restart:

#
# Define gamma as local (don't like doing this!)
#
  local gamma:

  A__h:= .18:        beta__1:= 0.8354e-1: psi:= 0.7258e-1:
  mu__r:= 0.51e-1:   sigma:= .165:        alpha:= .65:
  xi:= 0.5e-1:       gamma:= .131:        A__b:= .7241:
  beta__o:= 0.65e-1: mu__b = 0.17e-1:

  odes:=  diff(S(T), T) = A__h-psi*beta__1*S(T)*G(T)-mu__r*S(T),
          diff(G(T), T) = psi*beta__1*S(T)*G(T)-sigma*psi*beta__1*R(T)*B(T)-(alpha+xi+mu__r)*G(T),
          diff(A(T), T) = alpha*G(T)-(gamma+mu__r)*A(T),
          diff(R(T), T) = gamma*A(T)+sigma*psi*beta__1*R(T)*B(T)- mu__r*R(T),
          diff(C(T), T) = A__b - psi*beta__o*C(T)*I(T)- mu__b*C(T),
          diff(B(T), T) = psi*beta__o*C(T)*I(T)- mu__b*B(T):
  ics:= S(0)=100, G(0)=190, A(0)=45, R(0)=20, C(0)=35, B(0)=25:

#
# Solve system
#
  ans:= dsolve([ odes, ics ], numeric);

Warning, The use of global variables in numerical ODE problems is deprecated, and will be removed in a future release. Use the 'parameters' argument instead (see ?dsolve,numeric,parameters)

 

(1)

 

Download odeProb.mw

 

 

maybe?


You can adjust the decay rate by adding a parameter to the exponent in the definition of the function g(). I used 0.5, just for illustration.

restart;
g:=z->exp(-0.5*(z-24*floor(z/24)));
f:=(t)->g(t)*(t-24*floor(t/24))*150;
plot(f, 0..100);

proc (z) options operator, arrow; exp(-.5*z+12.0*floor((1/24)*z)) end proc

 

proc (t) options operator, arrow; 150*g(t)*(t-24*floor((1/24)*t)) end proc

 

 

 


 

Download staircase2.mw

It seems as if you want to add a "staircase" component to a function whcih already exists (I might be wrong!)

You might want to consider the "toy" example in the attached which adds 150 to the function 20*sin(z), at z=24, 48, etc

restart;
g:=z->20*sin(z);
f:=(t)->g(t)+(floor(t/24))*150;
plot(f, 0..100);

proc (z) options operator, arrow; 20*sin(z) end proc

 

proc (t) options operator, arrow; g(t)+150*floor((1/24)*t) end proc

 

 

 


 

Download staircase.mw

In the attached, I'm pretty sure that I am comparing the results using Physics:-Version(426) and Physics:-Version(429), although the message returned by the Physics:-Version() command is somewhat ambiguous.

In any case the two execution groups return different answers.

The first, using Physics:-Version(426) returns what I would consider to be the "correct" answer and the pdetest() command comfirms it (ie returns 0).

The second using Physics:-Version(429), returns the answer in your original post, which I agree is wrong!

Are you some kind of official beta-tester for Maplesoft??

restart;
Physics:-Version(426);
restart;
Physics:-Version();
pde:=diff(u(x,t),t$2)=4*diff(u(x,t),x$2);
ic:=u(x,0)=0,D[2](u)(x,0)=sin(x)^2;
bc:=u(-Pi,t)=0,u(Pi,t)=0;
sol:=pdsolve([pde,ic,bc],u(x,t));
pdetest(sol,pde);
simplify(diff(rhs(sol),t$2)-4*diff(rhs(sol),x$2));

Warning, this package updates content shipped in a standard Maple install.  Use the 'restart' command to clear your session before using these commands.

 

Kernel(`The "Physics Updates" version "426" is installed but is not active. The active version of Physics is within the library C:\Users\TomLeslie\maple\toolbox/2019/Physics Updates/lib\Physics Updates.maple, created 2019, September 23, 10:13 hours`), [`The "Physics Updates" version "426" is installed but is not active. The active version of Physics is within the library C:\Users\TomLeslie\maple\toolbox/2019/Physics Updates/lib\Physics Updates.maple, created 2019, September 23, 10:13 hours`]

 

`The "Physics Updates" version "426" is installed but is not active. The active version of Physics is within the library C:\Users\TomLeslie\maple\toolbox/2019/Physics Updates/lib\Physics Updates.maple, created 2019, September 23, 10:13 hours`

 

diff(diff(u(x, t), t), t) = 4*(diff(diff(u(x, t), x), x))

 

u(x, 0) = 0, (D[2](u))(x, 0) = sin(x)^2

 

u(-Pi, t) = 0, u(Pi, t) = 0

 

u(x, t) = Sum((1/2)*(Int(sin(n*x)*sin(x)^2, x = -Pi .. Pi))*sin(2*n*t)*sin(n*x)/((Int(sin(n*x)^2, x = -Pi .. Pi))*n), n = 1 .. infinity)

 

0

 

0

(1)

restart;
Physics:-Version(429);
restart;
Physics:-Version();
pde:=diff(u(x,t),t$2)=4*diff(u(x,t),x$2);
ic:=u(x,0)=0,D[2](u)(x,0)=sin(x)^2;
bc:=u(-Pi,t)=0,u(Pi,t)=0;
sol:=pdsolve([pde,ic,bc],u(x,t));
pdetest(sol,pde);
simplify(diff(rhs(sol),t$2)-4*diff(rhs(sol),x$2));

Warning, this package updates content shipped in a standard Maple install.  Use the 'restart' command to clear your session before using these commands.

 

Kernel(`The "Physics Updates" version "429" is installed but is not active. The active version of Physics is within the library C:\Users\TomLeslie\maple\toolbox/2019/Physics Updates/lib\Physics Updates.maple, created 2019, September 23, 10:14 hours`), [`The "Physics Updates" version "429" is installed but is not active. The active version of Physics is within the library C:\Users\TomLeslie\maple\toolbox/2019/Physics Updates/lib\Physics Updates.maple, created 2019, September 23, 10:14 hours`]

 

`The "Physics Updates" version "429" is installed but is not active. The active version of Physics is within the library C:\Users\TomLeslie\maple\toolbox/2019/Physics Updates/lib\Physics Updates.maple, created 2019, September 23, 10:14 hours`

 

diff(diff(u(x, t), t), t) = 4*(diff(diff(u(x, t), x), x))

 

u(x, 0) = 0, (D[2](u))(x, 0) = sin(x)^2

 

u(-Pi, t) = 0, u(Pi, t) = 0

 

u(x, t) = t*sin(x)^2

 

-8*t*(2*cos(x)^2-1)

 

-16*t*cos(x)^2+8*t

(2)

 

Download WaveEq.mw

just because alternatives are good!


 

restart;
with(inttrans):
alias( U__1(s)=laplace(u__1(t), t, s),
       U__2(s)=laplace(u__2(t), t, s),
       Y__1(s)=laplace(y__1(t), t, s),
       Y__2(s)=laplace(y__2(t), t, s)
     ):
sys:=[ 2*diff(y__1(t),t)=-2*y__1(t)-3*y__2(t)+2*u__1(t),
       diff(y__2(t),t)=4*y__1(t)-6*y__2(t)+2*u__1(t)+4*u__2(t)
     ];
lsys:=eval( laplace~(sys, t,s), [y__1(0)=0, y__2(0)=0]);

[2*(diff(y__1(t), t)) = -2*y__1(t)-3*y__2(t)+2*u__1(t), diff(y__2(t), t) = 4*y__1(t)-6*y__2(t)+2*u__1(t)+4*u__2(t)]

 

[2*s*Y__1(s) = -2*Y__1(s)-3*Y__2(s)+2*U__1(s), s*Y__2(s) = 4*Y__1(s)-6*Y__2(s)+2*U__1(s)+4*U__2(s)]

(1)

eq1:=isolate( subs(isolate(lsys[1], Y__2(s)),lsys[2]),Y__1(s));
eq2:=simplify(isolate( subs(eq1, lsys[2]),Y__2(s)));

Y__1(s) = (-2*U__1(s)+4*U__2(s)-(2/3)*s*U__1(s))/(-(2/3)*s^2-(14/3)*s-8)

 

Y__2(s) = ((2*s+6)*U__1(s)+4*U__2(s)*(s+1))/(s^2+7*s+12)

(2)

G_11:= simplify(eval(eq1, U__2(s)=0))/U__1(s);
G_12:= simplify(eval(eq1, U__1(s)=0))/U__2(s);
G_21:= simplify(eval(eq2, U__2(s)=0))/U__1(s);
G_22:= simplify(eval(eq2, U__1(s)=0))/U__2(s);

Y__1(s)/U__1(s) = 1/(s+4)

 

Y__1(s)/U__2(s) = -6/(s^2+7*s+12)

 

Y__2(s)/U__1(s) = 2/(s+4)

 

Y__2(s)/U__2(s) = 4*(s+1)/(s^2+7*s+12)

(3)

 


 

Download lap.mw

Isn't this a "straightforward" z-gradient scheme with a defined zrange, as in the attached (which has the additional benefit of working in Maple 2015).

Or am I missing something??

PS Plot colours render a lot better in Maple than they do on this site!
 

  restart:
  interface(version);
  with(plots):
  with(Statistics):
  randomize():
  N := 10:
  P := 3:
  A := Sample(Uniform(0, 1), [N, P]):
  C := CorrelationMatrix(A);
  matrixplot
  ( C,
    heights=histogram,
    axes=frame,
    gap=0.25,
    colorscheme=[ "zgradient",
                  [ "Blue", "White", "Red" ],
                  zrange=-1..1
                ],
    orientation=[0, 0, 0],
    lightmodel=none,
    tickmarks=[[seq(i+1/2=A||i, i=1..P)], [seq(i+1/2=A||i, i=1..P)], default],
    labels=[("")$3]
  );

`Standard Worksheet Interface, Maple 2015.2, Windows 7, December 21 2015 Build ID 1097895`

 

C := Matrix(3, 3, {(1, 1) = 1.0, (1, 2) = -.1563106801276543, (1, 3) = -0.6792568803622306e-1, (2, 1) = -.1563106801276543, (2, 2) = 1.0, (2, 3) = -.37067762598813236, (3, 1) = -0.6792568803622306e-1, (3, 2) = -.37067762598813236, (3, 3) = 1.0}, datatype = float[8])

 

 

 


 

Download zgrad.mw

First 8 9 10 11 12 13 14 Last Page 10 of 123