Carl Love

Carl Love

28035 Reputation

25 Badges

12 years, 321 days
Himself
Wayland, Massachusetts, United States
My name was formerly Carl Devore.

MaplePrimes Activity


These are replies submitted by Carl Love

@Kitonum Certainly there are some nontrivial combinations of the parameters Pr, lambda, and s for which Maple can "do" the integral that dsolve returns. One that works is Pr= 1, lambda= 1, s= 1.

@student_md You forgot about making X a Vector. Do something like X:= tau-> Vector(n, i-> xi[i](tau)).

I think that using x for gamma is a bad and confusing choice! How about gamma__0 (the 0 could be anything else)?

I don't think that there's hope for a symbolic solution because of the summation in f. I don't think that there'll be any problem getting a numeric solution (except perhaps if M is singular for some parameter values).

Although it's not an error, in the f definition you should change (i,j)-> to i->. Maple syntax allows unused extra parameters, but clearly in this situation there's no reason, plus it's confusing to the reader.

This is totally not your fault---just an unfortunate aspect of the various languages in use here (English, Maple jargon, and Maple code)---that you cause confusion both for yourself and others by using the word "table". There is perhaps no other word that's commonly used in all three of those languages that has as many different meanings.

Rather than listing the numerous possible meanings, perhaps it'd be easier if I asked you to simply define precisely what you mean by "table". And perhaps it'd be even better if you said precisely what you intend to do with a "table" once I show you how to obtain one from your graph.

@Rohith The following two situations are completely different, and require totally separate treatment:

  1. The function appears as a factor of the denominator, e.g., ln(x)/(a*sin(x))^2
  2. The function appears as a term of the denominator, e.g., 1/(a+2*sin(a)^2)

So, which of these is your situation?

@vv wrote:

  •  If you use the compiler, why do you keep the evalhf part?

The evalhf procedure and the compiled procedure do different things. The evalhf procedure is not used for the Binomial computation. It would be used, for example, for Mul(k^2, k= 1..10) mod 7. Its first argument is a dynamic procedure. I don't think that it's currently possible to pass that to compiled code.

  •  BTW, the compiler does not accept max with >2 arguments, you must change this.

Apparently this varies among compilers, because I had no trouble with it. Nonetheless, I removed it so that it will work for all.

  • without try/catch...

What??? What problem do you have with try/catch? Surely it must be used with any compiled code because it is so easy for it to error.

  • I also don't thing that recursion is useful (I have mentioned this in a previous comment).

And I mentioned in a comment below from about a week ago that I wasn't using recursion. My code does not now, nor did it ever, use recursion. I mentioned recursion before, but I was using the word incorrectly.

@jga I just posted new code and worksheet in the Reply three spaces above to address all the issues that you and VV brought up. Apparently my compiler is different than yours and VV's and mine accepts the three-argument max. Nonetheless, I removed it so that the code will work for all. (Reading ?Compiler,Compile, it certainly seems possible that there are various compilers in use.)

Regarding your time testing: For your time test of Maple's own binomial: Timing values that small are nearly all noise and have no relation to times that you'll experience in actual usage. They are smaller than the finest resolution of the "instrument" used to measure them. Don't believe the three decimal places. To test accurately, you need to repeat the operation many thousands of times using different input. I show an example below.

Regarding your time test of my Binomial: You are effectively only measuring the time required for the compiler to do its thing (returning an error in this case). The time for the Maple-level Binomial code is so small in comparison that it wouldn't change the last decimal place on those times.

Here's how to do realistic timings (this example is also in the worksheet):

(N,K,P):= 
   #3 64 x 64 matrices with roughly 4-digit entries
   rtable((1..2^6)$2, random(0..2^12)),
   rtable((1..2^6)$2, random(0..2^11)),
   nextprime~(rtable((1..2^6)$2, random(0..2^10)))
:
B1:= CodeTools:-Usage(
   Binomial~(N,K) mod~ P
):
memory used=16.22MiB, alloc change=8.00MiB, cpu time=109.00ms, real time=123.00ms, gc time=0ns

B2:= CodeTools:-Usage(
   binomial~(N,K) mod~ P
):
memory used=1.00GiB, alloc change=80.00MiB, cpu time=3.69s, real time=3.19s, gc time=1.33s


LinearAlgebra:-Equal(B1,B2);
              true

 

 

 

 

@vv The ability to diagnose the situation would be improved greatly by a positive case: that is, some input to isolve that returns with _SolutionsMayBeLost set. What it's set to doesn't matter much.

I based what I said largely on the fact that there have been a great many changes to solve since Maple 10, whereas the style of the isolve code suggests that it hasn't been updated in a very long time.

Do (and prepare for disappointment):

showstat(isolve);
showstat(`isolve/isolve`);
trace(solve):
sys:={x+y=10,x^2-y^2+z^2=1};
isolve(sys, n);

So, despite over a 100 lines of code, pretty much all it does is pass the problem to solve and then make a minor adjustment to what solve returns. I think that the setting of _SolutionsMayBeLost is entirely up to solve.

@digerdiga I didn't say that it was only true for 0 < t < Pi; I just said that it was true for those t, which left open the situation for other t. I was just providing an example of a single thing that you could change in your command that would make the simplification valid. Indeed, the interval could be stretched to -Pi < t <= Pi.

@vv You can replace simplify(convert(..., exp)) with combine(...).

@brian bovril Yes, your iteration is correct. But I don't want you to get the wrong impression that the reason to use x=y on the first iteration is because it allows the simplification of ln(1) = 0. That's just a side benefit. The reason is that x+ln(x) is nearly the identity function on your interval of interest, so x=y is already close to the solution.

Also, three iterations are required to be certain of getting the three decimals places that you asked for over the domain that you asked for. You show two iterations. Two iterations will get three decimal places at the lower end of the domain, but not the higher end.

@Rohith The commands selectremove, and selectremove only search the top-level operands of an expression. This job requires indets. I'll give more details later.

[In the following, W(x) (without boldface) refers to the mathematical Lambert W function rather than to my procedure W. On MaplePrimes, I always use boldface for inline references to literal Maple code.]

As VV has astutely observed, my first "fold" W(x,x) makes x the starting point for the iteration, and this convenient choice eliminates the ln from the first iteration. That convenience was not the reason that I chose this.

Since W(x) is the inverse of x*exp(x), it's fairly close to ln(x). Both ln(x) and W(x) are steep for 0 < x < 1. Steepness causes difficulty for Newton's method. But the function G(x) = W(exp(x)) is nearly the identity function---no steepness. It's inverse is x+ln(x) (as can be verified with solve), and thus the solution x of x+ln(x)=ln(y) is W(y). So, it's this equation rather than the more usual x*exp(x) = y that I solve with Newton's and/or Householder's method. Because of the closeness to the identity function, x0 = y is a great choice for the initial guess.

x + ln(x) = ln(y)  =>  x+ln(x/y) = 0. I define g:= (x,y)-> x+ln(x/y). The Householder method of order 1 is the ordinary Newton's method. The Householder method of any order is given by this simple procedure:

Householder:= proc(F, ord::posint)
local x, y;
   unapply(x+ord*D[1$ord-1](1/F)(x,y)/D[1$ord](1/F)(x,y), (x,y))
end proc:    

where F is a two-parameter procedure such that F(x,y) = 0 is the equation that we want to solve for x. Doing Householder(g, 3) puts a reasonably sized expression on my screen. Then I meticulously hand-optimized it by removing common subexpressions, subfactoring, etc. Why did I choose order 3? Because it seemed to have a nice balance between the simplicity of the iterative function and the number of iterations required for the desired accuracy. I haven't  checked, but it may minimize the total operation count.
 

 

@acer It's fine by me. Duplicate answers are unavoidable. You know what they say about great minds. I was surprised that your initial Answer didn't mention eval.

If a function appears in the expression without any coefficient, would you like its coefficient listed as 1, or would you rather that it not be listed at all?

First 294 295 296 297 298 299 300 Last Page 296 of 708