Carl Love

Carl Love

28055 Reputation

25 Badges

13 years, 1 days
Himself
Wayland, Massachusetts, United States
My name was formerly Carl Devore.

MaplePrimes Activity


These are replies submitted by Carl Love

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?

In your point 2, you have made an unsubstantiated claim of a serious error in Maple. You have had all day to post some evidence, but you have not. All claims of errors, even if unsubstantiated, cause damage to the reputation of the software and the company that makes it. As such, it is unethical to make those claims, just like it would be unethical to spread unsubstantiated claims that an associate of yours had a contagious disease.

So, I give you 12 more hours to post some evidence. If you won't I'll delete this thread, unless some other respected moderators can tell me some good reason not to.

@phbmc Yes, this is a 2D-Input problem that I've encountered before. It can't handle parameters whose type specifications refer to other parameters. So, if you make the typespec for simply L::list, it'll work in 2D. You may, if you wish, put an error check in the body of the procedure:

     if nops(L) < 2*c then error "List not long enough." fi;

Maple, oddly enough, doesn't care if you pass unneeded arguments to a procedure. It'll just ignore them (unless the code specifically looks for them). So, you can just pass in as a third argument (but I think that it's a weird thing to do). But, I don't see how it can be done without passing the list of items to be permuted. Does that make sense to you?

If you have a collection of related procedures, they can be (and probably should be) collected together into a module. That can easily depend on m, and I wouldn't even have any stylistic objections.

@zarara I know that it seems too easy to be true, but the limit command that VV shows does it all in this case, which is a super-easy problem. That is, if you're willing to accept the output of a computer program as "proof"...which I am. The big-O-ness is equivalent to that limit not being infinity, when that limit exists

By "snap" do you mean to click on, as with a mouse or other pointing device?

[This Reply is not related to binomial-mod-prime computations, and is intended only for expert-level module writers who write some numerical or semi-numerical code. It relates to a technique that I used in the module above.]

I just developed a coding technique---shown in the module above---that may be novel. It allows lexical variables (in this case module locals) to be used in evalhf'able and compilable procedures. The evalf'able and compilable procedures are defined containing a dummy variable _M. The values of _M are potentially machine dependent (at least theoretically), so they are computed in the ModuleLoad and subs into the procedures. The same thing could be done for any values that are not known until runtime.

I am saddened by the lack of enthusiasm for this Post. I believe that the code above is the only implementation in Maple of any exact[*1] statistical test of independence in a contingency table (the same property that is tested by the inexact Statistics:-ChiSquareIndependenceTest). Some modern statisticians recommend that given that we have the necessary mathematical software, that exact tests always be used for small (the sum of the entries in the table), say n < 1000. Furthermore the code above works for tables of any number of cells (i.e., the number of entries in the matrix) whereas some other programs that do Fisher's exact test are limited to two rows, two columns, or both.

[*1]"Exact" in this context refers to the fact that the computations are done with the relevant discrete distributions, the multinomial (generalized hypergeometric), rather than with a continuous distribution such as Chi-squared, which is only asymptotically equivalent to the discrete distributions as approaches infinity. "Exact" does not refer to the fact that the computations are done in exact rational-number arithmetic; however, the code above does that also. 

@vv I understand now why there was that factor-of-4 time difference between our codes in the last test that you showed: I was trying to do too much with just one variety of Mul, which handled arbitrary expressions akin to mul. I improved it by adding a separate Mul for ranges of integers, such as factorials, and I made it compiled for the appropriate cases. Please try it out.

I may have used the word "recursive" incorrectly earlier, leading to some confusion on your part. I simply meant that the procedure Binomial makes a single call back to itself to handle the smaller binomials. I didn't mean that the Lucas code was being called again.

First 295 296 297 298 299 300 301 Last Page 297 of 709