Allan Wittkopf

240 Reputation

7 Badges

18 years, 250 days

MaplePrimes Activity


These are replies submitted by Allan Wittkopf

I was unable to repeat the error you observed with the 'problem-1' worksheet in any of Maple 9, 9.5, 10, or the version currently in beta. I *did* however get a message when attempting to obtain the solution: Warning, computation is being performed near the boundary of the current precision, suggest increasing Digits to approximately 17 or efficiency may be degraded then the computation failed with: Error, (in dsolve/numeric/bvp) unable to achieve requested accuracy due to round-off error of approximately -.2e-2. Suggest increasing 'Digits' or selecting a higher 'abserr' tolerance I then bumped Digits up to 14 (use 15 if you are not using Windows), and dropped abserr to 1e-2, and was able to get a solution with no warnings. See Download 233_2847_problem-1.mws
View file details When I look at the solution, I see that some components have magnitudes of 1e6, so abserr being set to 1e-2 provides 8-9 digits of accuracy in the solution - if you need more, you will need to set Digits to a larger value (say 20), and give the computation a great deal more time, as it will need to work with non-hardware floating point. So questions: 1) What version of Maple are you using, and on what O/S? 2) Do you have a .mapleinit (non-Windows) or maple.ini (Windows) file that sets up your environment. If you do, what global variables do you define? I am looking at the second worksheet now. Sincerely, Allan
I was unable to repeat the error you observed with the 'problem-1' worksheet in any of Maple 9, 9.5, 10, or the version currently in beta. I *did* however get a message when attempting to obtain the solution: Warning, computation is being performed near the boundary of the current precision, suggest increasing Digits to approximately 17 or efficiency may be degraded then the computation failed with: Error, (in dsolve/numeric/bvp) unable to achieve requested accuracy due to round-off error of approximately -.2e-2. Suggest increasing 'Digits' or selecting a higher 'abserr' tolerance I then bumped Digits up to 14 (use 15 if you are not using Windows), and dropped abserr to 1e-2, and was able to get a solution with no warnings. See Download 233_2847_problem-1.mws
View file details When I look at the solution, I see that some components have magnitudes of 1e6, so abserr being set to 1e-2 provides 8-9 digits of accuracy in the solution - if you need more, you will need to set Digits to a larger value (say 20), and give the computation a great deal more time, as it will need to work with non-hardware floating point. So questions: 1) What version of Maple are you using, and on what O/S? 2) Do you have a .mapleinit (non-Windows) or maple.ini (Windows) file that sets up your environment. If you do, what global variables do you define? I am looking at the second worksheet now. Sincerely, Allan
When you just enter, for example, u_xy, it is by default not atomic. When you highlight it and make it atomic, then it becomes atomic. So if you do: u := cos: u_xy (Do not make this atomic) returns: cos_xy u_xy (highlight and make atomic) returns: u_xy So in the first case, the assignment, you assigned atomic u_xy to '2' meaning regular u_xy is unassigned. Then when you entered u := cos and u_xy, you did not make the u_xy atomic, so it returned cos_xy. If you had made it atomic before hitting return you would have received 2 instead. In the next bit, when you assign v_y=12 to 'test', and made v_y atomic, then the reference to v_y assigned in test is the atomic one. If instead you did the same thing as for the u example, just entered v_y without making it atomic, then you'd have received sin_y instead. Basically once an object is made atomic it becomes an immutable name that looks like what was made atomic. You then need to either assign it to a variable (at the time you made it atomic) for use, or make atomic any subsequent references that you want to also refer to the atomic object manually every time.
Yes, I believe that this didn't work quite right for Maple 10.00 and 10.01. If you upgrade to 10.04 or 10.06 I think that this will take care of the problem. - Allan
I just tried this in Maple 10.04, and had no problem (I am attaching the file). Could you let me know what version of Maple you are using? Just type 'kernelopts(version);' There were some problems in earlier point releases of Maple. Assuming you are using Maple 10.04 or later, there's one other thing. When you highlight the u_xy, then right click and go to "Arrange", what items other than "Atomic" do you see in the submenu? It should be "Subscript Up" and "Subscript Down". If not then it is possible that the selection may not be correct (alas you need to be very careful not to include any spaces to the left or right of the u_xy to get the right submenu). - Allan
Could you provide a bit more information on the problem you are solving? How does the main PDE depend on the solution of the other PDE? Is the other PDE also elliptic? - Allan
Maple does not currently have an elliptic PDE solver available.
Hello, The first thing I'd check is that your known function returns unevaluated if the input is not numeric. This is a requirement. For example: > f := unapply(x^2,x,numeric=x): > f(2); 4 > f(x); f(x) > eval(%,x=2); 4 Barring that, I could probably be more help if a worksheet demonstrating the problem could be posted. - Allan
Atomic objects are described in ?2DMathDetails, the last section. You can do: n := 1: Then later: f_n then highlight and make atomic and this will still be 'f sub n' not 'f sub 1'. The assignment to 'n' only affects table references. You have the same deal with f_3 - once made atomic both 'f' and '3' are locked. - Allan
By default, a subscript is a table reference, and f_max refers to the 'max' entry of the table 'f'. You can change this on an instance by instance basis by making the f_max atomic as follows: 1) Highlight the f_max 2) Right-click 3) Go to 'arrange' and select 'atomic' At this point, the f_max for which this was done is a subscripted object that is not a table reference, but can be thought of as the name 'f sub max'. - Allan
For starters, sorry for taking so long to respond. I was travelling and on vacation for most of July and some of August, and did not notice your response until recently. I'm still not entirely clear on what the issue is. The problem can be stated in terms of providing r(0), theta(0), r(10) and theta(10), or alternatively it can be written as an IVP in terms of r(0), r'(0), theta(0), theta'(0), but that is as many degrees of freedom as are available in the problem. Attempting to specify any more conditions will result in an overdetermined problem. You did mention that the 'time' is a flexible value, so it may be possible to solve this as an IVP, where you turn the thruster on until the conditions of being in the second orbit are satisfied, but for this, you need to construct a criteria to detect when you are in the desired second orbit. The issue here is that you are starting in the first orbit, and all you have available for control is the thruster (on/off) so you should not be able to have a completely free choice of time, and all parameters of the orbit, though it may be possible to get to an equivalent orbit in some unknown time (computed). If your thruster also had a direction associated with it, then it may be possible to get an exact orbit, but you'd need to vary both thrusts, and may still be unable to control at what time the satellite entered the orbit. I just uploaded a file, Download 233_OrbitIVP.mws , that shows that for this model, the provided conditions, and the choice of the thrust, it appears to only be possible to move the satellite into a higher but similar orbit. Perhaps the thrust model needs to be modified?
For starters, sorry for taking so long to respond. I was travelling and on vacation for most of July and some of August, and did not notice your response until recently. I'm still not entirely clear on what the issue is. The problem can be stated in terms of providing r(0), theta(0), r(10) and theta(10), or alternatively it can be written as an IVP in terms of r(0), r'(0), theta(0), theta'(0), but that is as many degrees of freedom as are available in the problem. Attempting to specify any more conditions will result in an overdetermined problem. You did mention that the 'time' is a flexible value, so it may be possible to solve this as an IVP, where you turn the thruster on until the conditions of being in the second orbit are satisfied, but for this, you need to construct a criteria to detect when you are in the desired second orbit. The issue here is that you are starting in the first orbit, and all you have available for control is the thruster (on/off) so you should not be able to have a completely free choice of time, and all parameters of the orbit, though it may be possible to get to an equivalent orbit in some unknown time (computed). If your thruster also had a direction associated with it, then it may be possible to get an exact orbit, but you'd need to vary both thrusts, and may still be unable to control at what time the satellite entered the orbit. I just uploaded a file, Download 233_OrbitIVP.mws , that shows that for this model, the provided conditions, and the choice of the thrust, it appears to only be possible to move the satellite into a higher but similar orbit. Perhaps the thrust model needs to be modified?
You need to use 'D' notation to specify the derivatives. Here is a workup of the solution for your problem with specific values chosen for the constants, and text plots for the solution at 3 time values:
> PDE := diff(u(x,t),t)=a*diff(u(x,t),x,x):
> IBC := {u(x,0)=T0, -k*D[1](u)(0,t)=q, -k*D[1](u)(d,t)=h*(u(d,t)-t1)}:
> vals := {a=1/10, T0=x^2, d=1, k=1, q=0, h=1, t1=0}:
> PDE := eval(PDE,vals);
                                         / 2         \
                       d                 |d          |
                PDE := -- u(x, t) = 1/10 |--- u(x, t)|
                       dt                |  2        |
                                         \dx         /

> IBC := eval(IBC,vals);
                     2
  IBC := {u(x, 0) = x , -D[1](u)(0, t) = 0, -D[1](u)(1, t) = u(1, t)}

> pds := pdsolve(PDE,IBC,numeric);
pds := module()
export plot, plot3d, animate, value, settings;
     ...
end module

> pds:-plot(t=0); 

 1 +                                                              AA  
   +                                                             A    
   +                                                           AA     
   +                                                         AA       
0.8+                                                       AAA        
   +                                                      AA          
   +                                                    AA            
   +                                                  AA              
0.6+                                                AA                
   +                                             AAA                  
   +                                            AA                    
   +                                         AAA                      
0.4+                                       AAA                        
   +                                    AAA                           
   +                                 AAA                              
   +                              AAA                                 
0.2+                          AAAA                                    
   +                      AAAA                                        
   +                 AAAAA                                            
   +          AAAAAAAA                                                
  -***********--+-+--+-+--+--+-+--+-+--+-+--+-+--+--+-+--+-+--+-+--+- 
               0.2          0.4         0.6          0.8           1  
> pds:-plot(t=1);

   +                                               AAAAAAAAA          
0.36                                           AAAA         AAAA      
   +                                        AAA                 AA    
0.34                                      AA                      AA  
   +                                    AA                            
0.32                                  AA                              
   +                                 A                                
0.3+                               AA                                 
   +                             AA                                   
0.28                           AA                                     
   +                         AA                                       
   +                        AA                                        
0.26                      AA                                          
   +                     A                                            
0.24                   AA                                             
   +                 AA                                               
0.22              AAA                                                 
   +            AA                                                    
0.2+         AAA                                                      
   +     AAAAA                                                        
0.18******-+-+--+-+--+-+--+--+-+--+-+--+-+--+-+--+--+-+--+-+--+-+--+- 
   0           0.2          0.4         0.6          0.8           1  
> pds:-plot(t=2);

   +                                 AAAAAAAAAA                       
0.27                              AAAA         AAA                    
   +                           AAA               AAA                  
   +                         AA                     AA                
   +                      AAA                        AA               
   +                    AA                             A              
0.26                  AA                                AA            
   +                AA                                   A            
   +             AAA                                      AA          
   +          AAA                                          A          
   +      AAAA                                              A         
0.25AAAAAA                                                   A        
   +                                                         AA       
   +                                                          A       
   +                                                           A      
   +                                                            A     
0.24                                                            A     
   +                                                             A    
   +                                                             AA   
   +                                                              A   
0.23--+-+--+-+--+-+--+-+--+--+-+--+-+--+-+--+-+--+--+-+--+-+--+-+--*- 
   0           0.2          0.4         0.6          0.8           1  

There is an example using derivative boundary conditions on the help page:
?pdsolve,numeric

- Allan
If your system is first order, then something as simple as the following will work:
> dsn := dsolve({diff(x(t),t)=y(t),diff(\                                      
> y(t),t)=-x(t),x(0)=0,y(0)=1}, [x(t),y(t)], numeric, output=operator);
dsn := [t = (proc(t)  ...  end proc), x = (proc(t)  ...  end proc),

    y = (proc(t)  ...  end proc)]

> xpr := eval(x(t)^2+y(t)^2,dsn[2..-1]);                                       
                                        2       2
                             xpr := x(t)  + y(t)

> eval(xpr,t=0);
                                      1.

> eval(xpr,t=1);
                                  1.000000253
For higher order systems, you can either re-write as a first order system, or substitute after converting the expression to 'D' notation:
> dsn := dsolve({diff(y(t),t,t)=-y(t),y(\                                      
> 0)=1,D(y)(0)=0}, [y(t)], numeric, output=operator);
dsn := [t = (proc(t)  ...  end proc), y = (proc(t)  ...  end proc),

    D(y) = (proc(t)  ...  end proc)]

> xpr := y(t)^2+diff(y(t),t)^2;
                                      2   /d      \2
                           xpr := y(t)  + |-- y(t)|
                                          \dt     /

> xpr := convert(xpr,D);
                                       2          2
                            xpr := y(t)  + D(y)(t)

> xpr := eval(xpr,dsn[2..-1]); 
                                       2          2
                            xpr := y(t)  + D(y)(t)

> eval(xpr,t=0);
                                      1.

> eval(xpr,t=1);
                                  1.000000253
Hope this helps, Allan
Sorry, the plots got scrambled by html conversion, trying again:

  +              AAAAAAA                               AAAAAA                 
0.35            AA      AA                           AA      AA               
  +           AA          AA                       AA          AA             
0.3          A             AA                     AA            AA            
  +         A                A                   A                A           
  +        A                  A                 AA                AA          
0.25      A                   AA               A                    A         
  +      AA                    AA             AA                    AA        
  +      A                      A             A                      A        
0.2     A                        A           A                        A       
  +    A                          A         AA                         A      
  +    A                          A         A                          A      
0.15  A                            A       A                            A     
  +  AA                            AA     AA                            AA    
0.1  A                              A     A                              A    
  + AA                              AA   AA                              AA   
  + A                                A   A                                A   
0.05A                                A   A                                A   
  +A                                  A A                                  A  
  *                                   AAA                                   A 
  *--+--+--+--+--+--+---+--+--+--+--+--*--+--+--+--+--+---+--+--+--+--+--+--* 
0               0.5              1             1.5              2             

   +                                                                  AAAAAA  
   +                                                              AAAAA       
   +                                                           AAA            
0.5+                                                         AA               
   +                                                      AAA                 
   +                                                    AAA                   
0.4+                                                  AA                      
   +                                               AAA                        
   +                                            AAA                           
   +                                        AAAA                              
0.3+                               AAAAAAAAA                                  
   +                           AAAA                                           
   +                        AAA                                               
0.2+                     AAA                                                  
   +                   AA                                                     
   +                AAA                                                       
   +              AAA                                                         
0.1+            AA                                                            
   +         AAA                                                              
   +    AAAAA                                                                 
  -******+--+--+--+--+--+--+--+--+--+--+--+--+--+--+--+--+--+--+--+--+--+--+- 
                 0.5             1             1.5             2              
1 2 3 4 5 6 Page 5 of 6