You appear to be using the quotient rule to compute the derivative of 4/x^3. When you did this the denominator should be v^2 = x^6, but you have x^5. This explains the difference between Maple's answer (-12/x^4) and yours (-12/x^3).

I would also point out that it is easier to differentiate expressions like 4/x^3 by writing the expression with a negative exponent: 4*x^(-3) and using the power rule: f''''(x) = 4*(-3)*x^(-3-1) = -12*x^(-4).

The work you included in your post made it easy to see your error. Thanks!

Doug

---------------------------------------------------------------------
Douglas B. Meade  <><
Math, USC, Columbia, SC 29208  E-mail: mailto:meade@math.sc.edu       
Phone:  (803) 777-6183         URL:    http://www.math.sc.ed

If you look at the end of the file you uploaded, it ends in the middle of a display.

It appears there are open Group and Presentation-Block tags. I closed these immediately before the final, incomplete, Ouput tag and removed the incomplete Output tag. I also closed the worksheet tag, but this is not necessary. The modified file opens in Maple 11 even though it still complains about being incomplete.

Here's the original ending of the file:

<Output>
<Text-field style="2D Output" layout="Maple Output"><Equation executable="false" style="2D Output" display="Pi... (but never terminated)

and here is what I changed it to:

</Output>
</Group></Presentation-Block></worksheet>

You can get the full file here:

View 178_4541_dxx-MOD.mw on MapleNet or Download 178_4541_dxx-MOD.mw
View file details

The student's output is huge. Do they really need to be doing exact arithmetic with more than 10000 digits?

I hope this is helpful to you, and to the student.

Doug

---------------------------------------------------------------------
Douglas B. Meade  <><
Math, USC, Columbia, SC 29208  E-mail: mailto:meade@math.sc.edu       
Phone:  (803) 777-6183         URL:    http://www.math.sc.ed

The best answer to your quesiton will depend somewhat on the details of the differential equations. I will illustrate some of the general ideas using some very simple examples. If there are not what you need, please give us some more details and we will try to find a workable solution.

I assume all of your differential equations are the same, just with a different parameter. (This example is so simple, this is NOT the way I would do THIS problem, but it illustrates the basic ideas to fulfill YOUR request. If your equations are simple, it would be easier to solve them in general.)

restart;
with( plots ):
ode := diff(x(t),t)=r*x(t);
                               d               
                              --- x(t) = r x(t)
                               dt              
ic := x(0)=1;
                                  x(0) = 1
R := [$1..100]/100:

s := rr -> dsolve( eval({ode,ic},r=rr), x(t), numeric );
s(1);
                    proc(x_rkf45)  ...  end;
s1(1)(1);
                    [t = 1., x(t) = 2.71828125045817615]

Now, to get the list of 100 points with the 100 different values of r, you can do the following:

p := seq( eval( [t,x(t)], s(r)(1) ), r=R );
          [1., 1.01005016708678386], [1., 1.02020134002697693], 

         ... (96 points omited)

            [1., 2.69123394153653804], [1., 2.71828125045817615]

These 100 points can now be plotted as follows:

plot( [p], style=point );

If you want to plot the 100 solutions on the interval [0,1], you could use the odeplot command as follows:

display( seq( odeplot(s(r), t=0..1 ), r=R ) );

I hope this will be helpful to you,

Doug

---------------------------------------------------------------------
Douglas B. Meade  <><
Math, USC, Columbia, SC 29208  E-mail: mailto:meade@math.sc.edu       
Phone:  (803) 777-6183         URL:    http://www.math.sc.ed

Excuse me, but didn't you ask about the same problem under the topic: Area Under Curves ?

Doug

---------------------------------------------------------------------
Douglas B. Meade  <><
Math, USC, Columbia, SC 29208  E-mail: mailto:meade@math.sc.edu       
Phone:  (803) 777-6183         URL:    http://www.math.sc.ed

My first thought is to point you to the select command. Here are some simple examples, if your situation is more complicated you will need to make some modifications to what I am about to demonstrate. If you need help, give us a better idea of your real need and we will see what we can do.

eq1 := x^2+x+1 = 0;
                                2            
                               x  + x + 1 = 0
s1 := solve( eq1, x );
                       1   1    (1/2)    1   1    (1/2)
                     - - + - I 3     , - - - - I 3     
                       2   2             2   2         
select( is, [s1], positive );
                                     []
eq2 := x^2+3*x-4 = 0;
                               2              
                              x  + 3 x - 4 = 0
s2 := solve( eq2, x );
                                    1, -4
s3 := fsolve( eq2, x );
                          -4.000000000, 1.000000000
select( is, [s2], positive )[];
                                      1
select( is, [s3], positive )[];
                                 1.000000000

If you want to do this in one step, you could define your own version of the solve command:

solvePOS := proc( eq, var )
  local s;
  s := solve( _passed );
  if s=NULL then
    s := fsolve( _passed );
  end if;
  select( is, [s], positive )[];
end proc:

solvePOS( eq1, x );

eq3 := expand( (x-1)*(x-2)*(x+1)*x ) = 0;
                           4      3    2          
                          x  - 2 x  - x  + 2 x = 0
solvePOS( eq3, x );
                                    1, 2

eq4 := x*exp(x) = 1.;
                                x exp(x) = 1.
solvePOS( eq4, x );
                                0.5671432904

I hope this is useful,

Doug

---------------------------------------------------------------------
Douglas B. Meade  <><
Math, USC, Columbia, SC 29208  E-mail: mailto:meade@math.sc.edu       
Phone:  (803) 777-6183         URL:    http://www.math.sc.ed

Robert's response is good, in particular the part about how 2D math represents the derivative. The way this has to be handled is not very attractive (at least to me). I much prefer the 1D approach, shown below:

ODE := diff(y(t),t)=x(t)+y(t)*(1-x(t)^2-y(t)^2);
                  d                     /        2       2\
                 --- y(t) = x(t) + y(t) \1 - x(t)  - y(t) /
                  dt                                       
SUBS := [ x(t)=r(t)*cos(Theta(t)), y(t)=r(t)*sin(Theta(t)) ];
           [x(t) = r(t) cos(Theta(t)), y(t) = r(t) sin(Theta(t))]

q := eval( ODE, SUBS );
/ d      \                                    / d          \                   
|--- r(t)| sin(Theta(t)) + r(t) cos(Theta(t)) |--- Theta(t)| = r(t) cos(Theta(t
\ dt     /                                    \ dt         /                   

                          /        2              2       2              2\
  )) + r(t) sin(Theta(t)) \1 - r(t)  cos(Theta(t))  - r(t)  sin(Theta(t)) /

Moreover, hitting this result with simplify yields some additional benefits:

simplify( q );
   / d      \                                    / d          \         /
   |--- r(t)| sin(Theta(t)) + r(t) cos(Theta(t)) |--- Theta(t)| = -r(t) \
   \ dt     /                                    \ dt         /          
                                        2              \
   -cos(Theta(t)) - sin(Theta(t)) + r(t)  sin(Theta(t))/

Either way, these are the basic tools to do what you want.

Doug

---------------------------------------------------------------------
Douglas B. Meade  <><
Math, USC, Columbia, SC 29208  E-mail: mailto:meade@math.sc.edu       
Phone:  (803) 777-6183         URL:    http://www.math.sc.ed

Do you want the four standard unit vectors in R4, or a general vector e[i]?

If it's the former, then you can use the following simple loop:

for i from 1 to 4 do
  e[i] := Vector(4,j->`if`(i=j,1,0));
end do;

I am not aware of any way to work with a standard unit vector indexed by the general i. You can probably define how this vector interacts with other vectors, but I'll not go there right now.

I hope this is helpful, or that others will show where I am wrong,

Doug

---------------------------------------------------------------------
Douglas B. Meade  <><
Math, USC, Columbia, SC 29208  E-mail: mailto:meade@math.sc.edu       
Phone:  (803) 777-6183         URL:    http://www.math.sc.edu/~meade/

restart; with( plots ):
p1 := plot3d( 2, y=-sqrt(4-x^2)..sqrt(4-x^2), x=-2..2, color=blue ):
p2 := plot3d( 4, y=-sqrt(16-x^2)..sqrt(16-x^2), x=-4..4, color=red ):
p3 := plot3d( sqrt(x^2+y^2), x=-5..5, y=-5..5, view=[DEFAULT,DEFAULT,2..4], color=green ):
display( p1,p2,p3, axes=boxed, transparency=0.7, style=patchnogrid );

You can get fancier, but I think this is straightforward and logical. Is it what you requested?

Doug

---------------------------------------------------------------------
Douglas B. Meade  <><
Math, USC, Columbia, SC 29208  E-mail: mailto:meade@math.sc.edu       
Phone:  (803) 777-6183         URL:    http://www.math.sc.edu/~meade/

Here is how I would approach this problem. First, s is the sum of the first n-1 integers:

s := n -> n*(n-1)/2;
                                     1          
                                n -> - n (n - 1)
                                     2

The list of ordered pairs (i,j) with i<j, i=1..n, and j=1..s can be constructed with

P := n -> [seq( seq( [r,c], c=r+1..s(n) ), r=1..n )];

P(3);
                          [[1, 2], [1, 3], [2, 3]]

P(4);
  [[1, 2], [1, 3], [1, 4], [1, 5], [1, 6], [2, 3], [2, 4], [2, 5], [2, 6], 
    [3, 4], [3, 5], [3, 6], [4, 5], [4, 6]]

To find the index in this list, use the member command with the optional third argument.

member( 1, [0,1,2], 'p' );
                                     true
p;
                                      2

These can all be put together in a single procedure, as follows:

LookUp := proc(e,S)
  local p;
  if not member(e,S,'p') then p:=0 end if;
  return p
end proc:

LookUp( [4,5], 5 );
                                      25
LookUp( [5,5], 5 );
                                      0

You might want to have a different behavior when the element is not in the set. I chose to have the LookUp procedure return 0.

Does this do what you have requested?

Doug

---------------------------------------------------------------------
Douglas B. Meade  <><
Math, USC, Columbia, SC 29208  E-mail: mailto:meade@math.sc.edu       
Phone:  (803) 777-6183         URL:    http://www.math.sc.edu/~meade/

Most of the functions available in the orthopoly package are now available in Maple without loading a package. In fact, the help for orthopoly states (at the bottom):

Note that this package will be made obsolete in a future version of Maple, and be replaced by a set of top-level functions such as HermiteH and LegendreP. See Also
ChebyshevT, ChebyshevU, GegenbauerC, HermiteH, JacobiP, LaguerreL, LegendreP, UsingPackages, with

You might want to check out these newer commands before you get too used to orthopoly.

Doug

---------------------------------------------------------------------
Douglas B. Meade  <><
Math, USC, Columbia, SC 29208  E-mail: mailto:meade@math.sc.edu       
Phone:  (803) 777-6183         URL:    http://www.math.sc.edu/~meade/

Yes, I do have a lot of experience with maplets. In my attempts to stay current with Maple's resources I have also learned a fair amount about embedded components.  I will post my comments in a separate blog entry in which I will also share with MaplePrimes a worksheet I created for some work I did to visualize sudoku puzzles. (This still needs a lot of work, and I am not opposed to accepting enhancements from others.)

For a lengthy response to acer's request, please see my blog entry titled Embedded Components, Maplets, and Sudoku Puzzles.

Doug
-------------------------------------------------------------------
Douglas B. Meade  <><
Math, USC, Columbia, SC 29208  E-mail: mailto:meade@math.sc.edu       
Phone:  (803) 777-6183         URL:    http://www.math.sc.edu/~meade/

There are many fewer embedded components (elements), and nothing like the Table feature in Maplets. Sure, you can mimic the effect using a table, but you have to put a separate TextArea component in each cell. This can be quite a nuisance. And, has recently been discussed in MaplePrimes, the same type of construction is possible with maplets.

Will the maplets package continue to evolve? Or, will future efforts be directed towards embedded components?

Doug

---------------------------------------------------------------------
Douglas B. Meade  <><
Math, USC, Columbia, SC 29208  E-mail: mailto:meade@math.sc.edu       
Phone:  (803) 777-6183         URL:    http://www.math.sc.edu/~meade/

I agree that it seems odd that it's not possible to change the values of a Table in a maplet, but the documentation certainly does indicate that the content of the table (including the header) is settable only initially.

An excerpt from the Maple help for Maplets,Elements,Table:

The Table element features can be modified by using options. To simplify specifying options in the Maplets package, certain options and contents can be set without using an equation. The following table lists elements, symbols, and types (in the left column) and the corresponding option or content (in the right column) to which inputs of this type are, by default, assigned. 
Elements, Symbols, or Types    Assumed Option or Content              
                                                                      
list of lists                  content of table body                  
list                           content of the table header            
range of positive integers     height and width options, respectively 
refID                          reference option                        
A Table element must contain exactly one TableHeader element and any number of TableRow elements. A Table element can contain a Font element to specify the font option. 
A Table element can be contained in a Maplet, BoxCell, or GridCell element, or Maplet element in a nested list representing a box layout. 
The following table describes the control and use of the Table element options. 
An x in the I column indicates that the option can be initialized, that is, specified in the calling sequence (element definition). 
An x in the R column indicates that the option is required in the calling sequence. 
An x in the G column indicates that the option can be read, that is, retrieved by using the Get tool. 
An x in the S column indicates that the option can be written, that is, set by using the SetOption element or the Set tool. 
Option     I  R  G  S  
                       
background x     x  x  
font       x        x  
foreground x     x  x  
height     x     x  x  
reference  x           
tooltip    x        x  
visible    x     x  x  
width      x     x  x  
                        

Can anyone explain this limitation?

Doug

---------------------------------------------------------------------
Douglas B. Meade  <><
Math, USC, Columbia, SC 29208  E-mail: mailto:meade@math.sc.edu       
Phone:  (803) 777-6183         URL:    http://www.math.sc.edu/~meade/

The use of loops in the initial implementation is the reason you see the noticeable increase in computational time as the length of the data increases. I will now show one way that avoids almost all loops.

We start with a clean slate, and define our data and the ranges that will be used for tabulating data.

restart;
input:=[144,144,108,94,92,93,96,96,96,96,96,99,95,101,101,97,98,93,95,100,96,94,95,98,97,95]:
R := [$0..255]:

The frequency table can be created in one line with

Gs := X -> [seq( nops(select(`=`,X,i)), i=0..255 )]/nops(X):
#Ga := X -> Array(0..255, i -> nops(select(`=`,X,i))/nops(X) );

The Gs procedure returns a list (with index 1 to 256). The Ga procedure does the same work, but returns an array (with index 0 to 255). Pick the one that you prefer; I have chosen to use Gs to more closely follow your approach. Note also that the Statistics package has a Tally command that could also be used, but the output is different enough that other changes would need to be made to construct the full frequency table.

TT := Statistics:-Tally( input );
 [144 = 2, 95 = 4, 94 = 2, 93 = 2, 92 = 1, 108 = 1, 96 = 6, 97 = 2, 98 = 2, 
   99 = 1, 101 = 2, 100 = 1]

GG := proc( input )
  local G;
  uses Statistics;
  G := Array( 1..256, 0 ):
  seq( assign('G'[lhs(e)]=rhs(e)/nops(input)), e=Tally(input) );
  return convert(G,list)
end proc:

If you do use Tally, it would probably be best to completely rework the rest of the encoding to take advantage of the specific output produced by Tally.

The real bottleneck in the encoding is the computation of the cumulative frequency table. Instead of trying to find a nifty implementation of this using recursion, etc., I will simply use the CumulativeSum procedure in the Statistics package. To be completely consistent with the original implementation, it is necessary to convert this table to rational numbers; if you are happy with floating point results, omit the conversion to rational in the definition of V.M

Encoding:=proc(input,R)
  local N,i,s,L,H,x,G,V;
  uses Statistics;
  N := nops(input);
  G := Gs(input);
  V := CumulativeSum( G );
  V := map( convert, V, rational );  # omit if floating point is OK
  s:=1: L:=0: H:=1:
  for x in input do
    if member(x,R,'i') then
      L:=L+s*V[i-1];
      s:=s*G[i];
    end if;
  end do;
  H:=L+s;
  return [L,H]
end proc:

Encoding( input, R );
     [1466769803662367936038630151511  1466769803662367936038630163175]
     [-------------------------------, -------------------------------]
     [1467733283092297866534393856144  1467733283092297866534393856144]

With my timings, your approach takes 0.1 CPU second and mine takes 0.04 CPU seconds. I did some tests with larger datasets, but it was not clear how much faster my code is.

Hoping this is useful,

Doug

---------------------------------------------------------------------
Douglas B. Meade  <><
Math, USC, Columbia, SC 29208  E-mail: mailto:meade@math.sc.edu       
Phone:  (803) 777-6183         URL:    http://www.math.sc.edu/~meade/

Did you take a look at the arrows and dirgrid options?

arrows affects the type of symbol drawn in the direction field. I'm guessing you might prefer arrows=line or arrows=mediumfill or arrows=curve.

dirgrid tells Maple the number of points, in each direction, to use to create the direction field. Increasing the number of points can give the illusion of solution curves, in some cases.

Here is my variation on your original code:

with( plots ):
with( DEtools ):
c:=-4:
FM1:=DEplot(diff(y(x),x)=-c/x^2,y(x),x=-10..10,y=-10..10,color=red,arrows=curve,dirgrid=[30,30]):
FM2:=DEplot(diff(y(x),x)=x^2/c,y(x),x=-10..10,y=-10..10,color=green,arrows=curve,dirgrid=[30,30]):
display(FM1,FM2);

Hoping this is helpful,

Doug

---------------------------------------------------------------------
Douglas B. Meade  <><
Math, USC, Columbia, SC 29208  E-mail: mailto:meade@math.sc.edu       
Phone:  (803) 777-6183         URL:    http://www.math.sc.edu/~meade/