Allan Wittkopf

240 Reputation

7 Badges

18 years, 283 days

MaplePrimes Activity


These are replies submitted by Allan Wittkopf

I am unable to repeat the Matlab result, and if I understand the problem correctly, for the given input system, the Matlab result must be wrong. Once again, back to the ODE argument. For the Matlab solution, it appears as though the value P(0) is in the region of 0.135. Now when I apply a shooting method with P(0) ~= 0.135, and D(P)(0)=0, the resulting plot is actually oscillitory, it does not vary smoothly as the Matlab plot does. In fact, if I obtain plots for P(0)=0.134 and P(0)=0.136, there is very little difference between them. Both functions are oscillitory, have 4 peaks and 5 valleys, and neither has D(P)(8e-8) even close to 0. Now since the ODE does not become undefined in the region z=0..8e-8, there is a unique solution for a given value of D(P)(0) and P(0) for this problem, and since there is little variation around P(0)=0.135, then the Matlab solution is not a solution of the IVP taken from z=0, therefore it is not the solution to the problem. The fact that there is such a dramatic difference makes me believe that either: 1) The solution in Matlab is an artificial solution introduced by the discretization. or 2) There is a very small error (perhaps a sign error) in the setup of the problem. Sorry I was unable to help, Allan
Just to confirm then, the *exact* problem you are solving is: A:=abs(3.0164e14*z-2.8071e7)*(-0.82500e-2*z+3e-9)^2: alpha_z:=3.0164e14*z-2.8071e7: beta_z:=9.9000e15*z-5.400e8: gamma_z:=6.6e9: # (taken as constant for now) bc := D(P)(0)=0, D(P)(8e-8)=0; deq := A*diff(P(z),z,z)-alpha_z*P(z)-beta_z*P(z)^3-gamma_z*P(z)^5=0; Is every coefficient and term in the above formulation correct to the specified number of decimal points? It takes some time to thoroughly analyze these problems, so getting it right would be helpful. Thanks, Allan
I don't seem to be able to upload the worksheet to the site, so here are the commands I used to demonstrate that the problem is single valued: # Define problem A:=abs(3.0164e14*z-2.8071e7)*(-0.82500e-2*z+3e-9)^2: alpha_z:=3.0164e14*z-2.8071e7: beta_z:=9.9000e15*z-5.400e8: gamma_z:=6.6e9: # (taken as constant for now) bc := D(P)(0)=0, D(P)(8e-8)=0; deq := A*diff(P(z),z,z)-alpha_z*P(z)-beta_z*P(z)^3-gamma_z*P(z)^6=0; # Generate a solution for D(P)(0)=0 and some arbitrary value for P(0) Digits := 30: dsn := dsolve({deq,D(P)(0)=0,P(0)=0},numeric,abserr=1e-20,relerr=1e-20): # Create a function that uses dsolve/numeric to compute the value D(P)(8e-8) fr> fcn := proc(P0) local val; if type(P0,'numeric') then dsn(initial=[0,P0/1e16,0]); try val := op([3,2],dsn(8e-8)); catch: val := Float(undefined); end try; val; else 'procname'(P0); end if; end proc: # Plot this function over several ranges: plot(1e13*fcn(z/1e6),z=-1..1); plot(1e10*fcn(z/1e3),z=-1..1); fcn(1e-3); fcn(1e-2); fcn(1e-1); fcn(1e0); fcn(1e1); # So from the plots, and from the above trend, it appears as though # the only value of P(0) that gives D(P)(8e-8)=0 is 0
The reason I did this is as follows: The ODE in question is leading linear, and the coefficients are purely polynomial in the interior of the region 0..8e-8. Since this is true, if the coefficient of the leading derivative of the problem is non-vanishing, then the initial value problem over the specified range is well defined. This suggests that use of a shooting approach will be a possibility for obtaining a solution. When I attempted to apply a shooting approach, varying the value of P(0), (with the incorrrect value for 'A'), I obtained a problem where the value of D(P)(8e-8) varied monotonically w.r.t. P(0), and the only value for which I could get D(P)(8e-8)=0 was for P(0)=0. Now that I have the corrected coefficient, the problem is behaving differently. I'll keep looking. - Allan
As it is defined: A:=abs(3.0164*10^14*z-2.8071*10^7)*(-0.82500e-2*z+3.*10^(9))^2; alpha_z:=3.0164*10^14*z-2.8071*10^7; beta_z:=9.9000*10^15*z-5.400*10^8; gamma_z:=6.6e9; # (taken as constant for now) bcs := D(P)(0)=0,D(P)(8e-8)=0; deq := A*diff(P(z),z,z)-alpha_z*P(z)-beta_z*P(z)^3-gamma_z*P(z)^6=0; the solution (unless I am mistaken) is single-valued over the speccified interval. If you look at the coefficient of diff(P(z),z,z) with the above definitions, you get: > lc := abs(.3016400000e15*z-28071000.00)*(-.82500e-2*z+3000000000.)^2; 15 8 lc := | 0.3016400000 10 z - 0.2807100000 10 | 10 2 (-0.0082500 z + 0.3000000000 10 ) > fsolve(lc=0,z=0..8e-8); 15 8 fsolve(| 0.3016400000 10 z - 0.2807100000 10 | 10 2 -7 (-0.0082500 z + 0.3000000000 10 ) = 0, z, 0 .. 0.8 10 ) So as far as I can tell, there should only be one solution, and that solution should be zero. Did I make the typo correction incorrectly in the provided value for A? - Allan
> Dear Allan! > > thanks for your post. I can imagine what amount of work that is, > and I am not sure if I want to endure this, and if I really could > do it. I know by now, that is essentially the java and the > dynamical loader which cause the problems and I was hoping for a > sort of a direct fix. Well Maple 6-7 had nothing to do with java at all, while Maple 8 only used java for Maplets, so really the problem is dynamical loading and the changes made on Linux to gcc. > I wonder what an 'unlimited' (in time) license does mean. Sure > enough, if one would have the sources then another make, make > install would do the trick. I see the legitimate interest and > need of maple soft to sell new products, but at some times I > would just be happy to be able to get an old version running > _painlessly_ It does run painlessly on an older O/S :-) Seriously, there's just no way to tell what will happen next on Linux. When the change to the behavior of glibc happened, I expect the developers worked to change usage of the relevant bits of glibc so that the binary would become independent of that change, but this could not be done until it was detected (which of course was after it happened) and corresponds to after Maple 8 was released. > Even maple 10 does not run out of the box on a decent linux > machine, like SuSE 10.1, and there is no Maple 11 around yet > which the should run (??) This is more glibc changes (sigh...) combined with binary incompatibilities wrt. those changes for java. Did you see the exceptionally good post about modifying the installer to comment out the LD_ASSUME_KERNEL in the installation script? > I think that I will possibly just run an old retired computer > without internet connection as a sort of old maple resort home > thing for testing programs across platforms. That's definitely the easier way (easier than the way I do it anyway...) - Allan
Another method of determining the point at which the solution attains a certain value is to use stop conditions: > dEq := diff(theta(t), t, t) = sin(theta(t)): > soln := dsolve({dEq, theta(0) = 0, (D(theta))(0) = 1}, theta(t),numeric, > stop_cond=[theta(t)-2]); soln := proc(x_rkf45) ... end proc > soln(10); Warning, cannot evaluate the solution further right of 1.4849919, stop condition #1 violated [t = 1.48499193139580199, theta(t) = 2.00000000000000044, d -- theta(t) = 1.95763036755598319] dt - Allan
Another method of determining the point at which the solution attains a certain value is to use stop conditions: > dEq := diff(theta(t), t, t) = sin(theta(t)): > soln := dsolve({dEq, theta(0) = 0, (D(theta))(0) = 1}, theta(t),numeric, > stop_cond=[theta(t)-2]); soln := proc(x_rkf45) ... end proc > soln(10); Warning, cannot evaluate the solution further right of 1.4849919, stop condition #1 violated [t = 1.48499193139580199, theta(t) = 2.00000000000000044, d -- theta(t) = 1.95763036755598319] dt - Allan
I have a set up that I use to do this, and I can currently run Maple V.4 through Maple 10 within this setup, but the process is somewhat lengthy. First the reason that this is a problem is not in fact the Linux kernel, but rather changes that were made to glibc to make Linux thread safe. These changes introduced a compatibility problem with any binaries that use 'errno'. So given that, it is not possible to run many older binaries on the more current distributions. One option for slightly older Linux distributions is to set: LD_ASSUME_KERNEL=2.2.5 and this causes libc to behave in a bit of a backward compatibility mode, but support for this is being phased out, and is not present in most newer distributions (FC2 or newer, Suse 9.? or newer, etc.) You can try to get a copy of the older libc and libm, then add to your LD_LIBRARY_PATH, but the dynamic linker is one thing you have little control over, and it will typically fail when you try to use an older version of libc. So in reality you really need to be running a Linux system with a libc version 2.2.5 or older to run many of these older applications. The only general alternative I have found is to construct a "chroot jail", using the following steps (the *very* short version of the instructions): 1) Install a nearly complete older Linux distribution in a directory on your machine (I use /usr/RH5.2 as I am using a RH5.2 distro for this purpose). The only safe way to do this is to copy it off of an existing functional distribution (you can't run an install from chroot). 2) Remove directories from the distro for which you actually want to use the directories on the core machine - this would include things like /dev and /proc 3) Add in bind mounts to remount existing mounted directories to mount points inside the /usr/RH5.2 subdirectories (you need to create the directories for the mount points). Examples: /bin/mount --bind / /usr/RH5.2/native /bin/mount --bind /usr /usr/RH5.2/nativeusr /bin/mount --bind /dev /usr/RH5.2/dev /bin/mount --bind /proc /usr/RH5.2/proc 4) chroot /usr/RH5.2 as root and very carefully replace the bits in /etc (which is really /usr/RH5.2/etc) with links to /native/etc (which is really /etc). This must be done from within the chroot, as otherwise the links will not work correctly). This would include /etc/passwd and the like - things you only want to maintain on the main machine. 5) su - At this point it will be as though you are running a RH5.2 system, using RH5.2 libc, libm, dynamic linker, X libraries, etc. Just run your programs from here (for 'X' do not forget to set your DISPLAY variable). Note that not everything works perfectly, for example you cannot open up an xterm from your chroot setup, because the ptys are not necc. compatible across that many kernel versions, but xmaple works just fine. - Allan
This bug is now fixed in Maple 11 (now in beta) As an aside, I have noticed that in your worksheets you often use a *very* fine initial mesh, even when it is not necessary. The resulting solution is obtained quite slowly. Use of a coarser initial mesh will provide a solution with the same accuracy much more rapidly in the cases where a solution can be obtained. In the example you attached to this post you use 'initmesh=8*128', and the computation runs 5 seconds on my machine. Instead choosing 'initmesh=128' produces a solution with the same accuracy in 0.7 sec, and a solution 1000 times more accurate in 0.73 sec. For a smaller problem (as in this example), the difference is not huge, but for larger problems with many equations, the difference can be much more substantial. Sincerely, Allan Wittkopf
This bug is now fixed in Maple 11 (now in beta) As an aside, I have noticed that in your worksheets you often use a *very* fine initial mesh, even when it is not necessary. The resulting solution is obtained quite slowly. Use of a coarser initial mesh will provide a solution with the same accuracy much more rapidly in the cases where a solution can be obtained. In the example you attached to this post you use 'initmesh=8*128', and the computation runs 5 seconds on my machine. Instead choosing 'initmesh=128' produces a solution with the same accuracy in 0.7 sec, and a solution 1000 times more accurate in 0.73 sec. For a smaller problem (as in this example), the difference is not huge, but for larger problems with many equations, the difference can be much more substantial. Sincerely, Allan Wittkopf
OK, this is a bug. Given where it occurs in the code, it will only arise when: 1) Continuation is in use 2) The problem requires mesh adaptation to reduce the error 3) No other steps occur between the completion of the continuation and the start of the mesh adaptation. The reason that this works when you set tolerance to 1e-2 is because mesh adaptation is not required to reach that error tolerance. This also shows up on Linux with a tolerance of 1e-5 for the same problem (Linux gets one more Digit of accuracy than Windows). Strangely enough, even dodging the bug (in the code), you still do not get a successful computation, as the requested tolerance is too small for hardware Digits (the solution is O(1e6)). This will be fixed for Maple 11. - Allan
OK, this is a bug. Given where it occurs in the code, it will only arise when: 1) Continuation is in use 2) The problem requires mesh adaptation to reduce the error 3) No other steps occur between the completion of the continuation and the start of the mesh adaptation. The reason that this works when you set tolerance to 1e-2 is because mesh adaptation is not required to reach that error tolerance. This also shows up on Linux with a tolerance of 1e-5 for the same problem (Linux gets one more Digit of accuracy than Windows). Strangely enough, even dodging the bug (in the code), you still do not get a successful computation, as the requested tolerance is too small for hardware Digits (the solution is O(1e6)). This will be fixed for Maple 11. - Allan
The problem with the first worksheet is a bug, and I notice that it is still in Maple 10 if you keep the tolerance at 1e-4. I will investigate.
The problem with the first worksheet is a bug, and I notice that it is still in Maple 10 if you keep the tolerance at 1e-4. I will investigate.
1 2 3 4 5 6 Page 4 of 6