Carl Love

Carl Love

28020 Reputation

25 Badges

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

MaplePrimes Activity


These are answers submitted by Carl Love

The conversion of an expression sequence to a list and vice versa is a triviality both in terms of coding and efficiency. Kitonum has shown the coding syntax. Regarding the efficiency, both the time and memory needed for the conversions are very small constants; they do not depend on the length of the data.

It is often better to stick with lists because of the indexing problems that can occur when an expression sequence has 0 or 1 elements.

This should not have been a separate Question, but since Tom already put an Answer here, I won't delete it.

To me, the input diff(x-1 = 0, x) is not GI (garbage input) nor is 1=0 GO (garbage output). As I was trying to explain in the other thread, everything makes sense if you take into account the quantifications that are implicitly present.

Maple allows for the simultaneous differentiation of both sides of an equation, and this is often useful. However, you could argue that this shouldn't be allowed (because of this very problem). And thus you could argue that diff(a = b, x) should be considered GI and an error for any a and b. If this restriction were to be imposed, I wouldn't object. It's trivial to use overload to impose that restriction:

restart:
diff_orig:= eval(diff):
unprotect(diff):
diff:= overload([
    proc(e::`=`)
    option overload; 
        error "Differentiation of equations not allowed."
    end proc, 
    diff_orig
]):
protect(diff, diff_orig):

 

The output 1=0 could be viewed as an algebraic equation or as a boolean relation. If it's an equation, then it's an equation whose solution set is empty. That's not garbage. If it's a relation, then it evaluates to falsefalse is not garbage.

The input and output considered together constitute a proof-by-contradiction that there does not exist any differentiable function f such that f(x) = x - 1 and f(x) = 0 (for all x in whatever domain of differentiation you want to consider). That's knowledge (albeit trivial), and knowledge isn't garbage.

You could do it in the way that you suggest if you want, but all the standard differential operators are in the VectorCalculus package: Gradient, JacobianHessian, and several others.

As a workaround, you can do this:

with(Physics):
Physics:-Setup(assumingusesAssume= false):
csgn(a - b[2]) assuming a > b[2];

The simplify and all the other assumptions are superfluous. You can put them if you want, but they make no difference.

 

The number of combinations that you're trying to generate is so large that they could not be listed by any computer (or any other method), ever. And it's not possible that there ever will exist a computer (or any other method) of doing this. We can count the number of such combinations like this:

Number(Iterator:-MultiCombination([2^16 $ 2^8], 2^16));
7820701318841130906613553459164687216584467046464667050034224507
  41484858519585179827709391651098624501404918769795637653274401
  57029711677787977640862663017116097938315215722518761350936266
  38516727105642689147056759932239793313894413127816886156554817
  76804366751588302164524267450481108817377756113734373353127476
  60436555681987839580644837879373377115558668241858194605420318
  57090864326621621084868745055201958236210246318706997513184396
  62089540815829177821487234240936427010441020860670200642172134
  59286677922608980438360909995794377076419365877451789293653200
  24753539064997228881393914518695160243181091267722013448011707
  00847339475028997519306136404337673374805865549565962012504690
  7222757907787394934267401909783755764225

evalf(%);
                                     723
                       7.820701319 10   

So, if every cubic meter in the entire Universe contained a computer, and they each generated one combination each nanosecond over the entire age of the Universe, then the progress they would've made on the job would still be far less than the progress made on the job of counting every atom in the Universe by counting a single atom.

So, you need to say what you want to accomplish by generating these and figure out a way to do it without generating them.

The way that you're using indexing (numbers appearing in square brackets appended to expressions) makes no sense to me, and the error message implies that it makes no sense to pdsolve either. I have no idea whether pdsolve will ultimately be able to solve your problem, but at this point it doesn't even recognize the problem as a PDE.

Suggestions:

  1. Don't use 2D input.
  2. Make your primary goal to get pdsolve to recognize your system as a PDE without error. Having the system prettyprinted on screen in a certain way should only be a distant secondary goal.
  3. Don't use PDEtools:-declare because its purpose is primarily to address cosmetic rather than syntactic issues.
  4. If you want 2 us (or any other variable), just call them u0 and u1
  5. Stop saying that you want to solve it by the perturbation method (or any other method). If you can't even input the problem, then the method chosen to solve it is at the moment irrelevant.

The cosmetic issues can be addressed afterwards. They are of course irrelevant if it can't even recognize what you've entered as a PDE.

By design, the command subs does not do any type of mathematical manipulation. The things appearing on the left sides of the substitutions must explicitly appear in the expression being changed. Understanding the mathematical relationship between exp(a) and exp(-a) is something that subs does not do (by design).

"Convert floats to rationals" means to convert any numbers with decimal points (aka "floats") into exact integers or fractions expressed as a ratio of integers (aka "rationals"). The command convert(..., rational) will do this. You can apply this command to the entire system of equations at once; there's no need to treat each number individually.

You are using the word "order" where I think you mean degree. Yes, it's true that your equations (after being reduced to polynomials) are degree 3. In the univariate case, that would mean 3 solutions per equation (counting possible multiplicities and solutions that must be discarded because they lead to 0 denominators). The multivariate case is more difficult to count. Without regard to your specifc equations, all that we can say is that the maximum number of solutions of a system of 2 degree-3 equations is 3*3 = 9.

Here's how to convert your system to rationals:

eq_VS_m4:= 132.790562 = -0.000588*VRi + 0.995381*VRr - 
    (1692.955600*VRi - 1747.964480*VRr)/(VRi^2 + VRr^2):
eq_VS_m5:= 0 = 0.000588*VRr + 0.995381*VRi + 
    (1747.964480*VRi + 1692.955600*VRr)/(VRi^2 + VRr^2):
Sys:= convert({eq_VS_m||(4..5)}, rational, exact);

The exact is optional.

And here is how the system can be reduced to a single quadratic:

Se:= eliminate(Sys, VRr);
      [ /        
Se := [{ VRr = - 
      [ \        

                               1                                      
  ------------------------------------------------------------ (4 VRi 
  33953752813744399363729225 VRi + 435026389377170333527426880        

  (8764251389865934151895845 VRi + 112422402706002345584760664))

  \    /                                2
   }, { 4367704732338974625213889705 VRi 
  /    \                                 

   + 111644550712064136430569347176 VRi

                                   \ ]
   + 713443278578559346984980083200 }]
                                   / ]

I chose to "eliminate" the variable VRrIn this particular case, I could just as well have eliminated VRi; however, it's not always that easy to switch. (Indeed, deciding on the order with which to eliminate variables is a major computational issue in computer algebra.) The result of eliminate is a list of two sets. The first set is equation(s) (just one equation here) solved for the variable(s) that you wanted to eliminate. The second set is expression(s) (assumed to be equated to 0) that remain to be solved. In this case, we see that the remaining equation is quadratic.

The Answers by Kitonum and Oliver are correct. I'd just like to provide an Answer with more theoretical context. And I realize that you realize that what you posted is more-or-less an anomaly rather than a bug. There is however a very broad CAS design issue that could be addressed.

I'll assume that you're at least somewhat familiar with the two main[*1] quantifiers of formal logic (although you may not have realized that these concepts had such formality). Their formal names are "the universal quantifier" and "the existential quantifier". There are a vast number of ways that these can be phrased in common language; the most-generic phrasings are "for every" or "for all" for the universal and "there exists" for the existential. In symbols, the universal usually appears as an upside-down uppercase A, and the existential as a mirror-imaged uppercase E. Quantifiers are used to bind[*2] variables to propositions (statements whose truth can be evaluated once their variables, if any, have been evaluated).

The problem with Maple (and it's also a problem in many CASs) is that variables are too often allowed to be free[*2] when they should be bound[*2].

Now let's consider your situation starting with your equation x+y = 1 in some different contexts (i.e., with some different variable bindings (the plural gerund of bind)). 

  1. In its most-usual algebraic context, this equation could be more formally stated as "There exists a pair (x,y) such that x+y = 1" and even more formally as "There exists a pair (x,y) in A x B such that x + y = 1" (where A and B are some specified sets such that the proposition is sensical and A x B means the Cartesian product).
  2. ​​​​​​In the context that you use the equation, it could be more formally stated "There exists a y in some domain D over which differentiation is defined such that for all x in D, x+y = 1."
  3. (Oliver's context) "There exists a function (or multi-function) y defined over some domain D over which differentiation is defined such that for all x in D, x+y(x) = 1."

So, these are mathematically different statements. In the context that you used, (2), obtaining 1 = 0 is not a paradoxical contradiction. Rather, it's a proof-by-contradiction contradiction that proves that there does not exist such a y.

[*1] Each of these quantifiers can be expressed in terms of the other.

[*2] The past participle of the verb bind is bound (e.g., "Variable x is bound in the expression int(f(x), x= a..b)."). The opposite of this usage of bound is the adjective free. This usage of bound mustn't be confused with the noun bound and its associated past participle bounded (e.g., "The upper bound of x in int(f(x), x= a..b) is b; x is bounded by a and b."). This can be quite tricky because, as the example sentences show, the two distinct mathematical meanings of bound can appear in identical contexts.

By default, Maple's Dirac has no assigned value at 0. You can assign the value in your preliminary code as a "remember table assignment":

Dirac(0):= 1

The old way that you mentioned still works in general, but your examples int and type are just unlucky choices. They are both special cases, but in completely diffferent ways.

Regarding type: It is a builtin command, meaning that it is not written in Maple and is "built in to" the kernel rather than stored in the Library. There is no way to list the code of builtin commands, and there never has been a way.

Regarding int: It has been recoded with a module wrapper. Your print method will not work with modules. But there is still an underlying procedure inside that module, and you can view it with showstat(int) or

interface(verboseproc= 3);
kernelopts(opaquemodules= false);
print(int:-ModuleApply);

When a command cmd has a module wrapper, the module's name is cmd and the procedure's name is cmd:-ModuleApply. When you're trying to read code, you usually want to set opaquemodules to false, as above, which gives module locals the same accessibility as exports (via the : - operator) except that the locals won't be affected by the with command. 

The codegen package is quite convenient for dynamic procedure generation, but unfortunately it is much older than keyword parameters and also (I think) default-value parameters. Here are two other methods for dynamic procedure creation: the subs method and the FromInert method. The subs method is sufficient for the vast majority of practical cases. The FromInert method can be unwieldy, but it can handle any case at all.

Method 1: The (undocumented) subs method: Although this usage of subs is undocumented, it has been in use at least 20 years and is used in Library code, so I don't expect it to change. We start with a "template" procedure, and create the dynamic procedure by making substitutions into it by using the syntax of the ordinary subs command (although the semantics are substantially different in undocumented ways).

Here is a template for a procedure similar to the one you showed. For brevity, I've omitted your parameters a and b.

#Procedure template:
P:= proc(  
    _rp1:= _rp1d, 
    _op2::_op2t:= _op2d, 
    {_kp3:= _kp3d, _kp4::_kp4t:= _kp4d}
) 
    (_rp1, _op2, _kp3, _kp4) #return value
end proc
:

Everything that will be changed by subs I've given a name beginning with underscore. It's not necessary to use an underscore, but they must be implicitly global symbols without assigned values. By implicitly, I mean that you can't force them to be global by using a :- prefix. Nor can you use unevaluation quotes '...' to ignore values that they've possibly been assigned.

Now we assign the desired parameter names, types, default values, and keywords:

#Create dynamic run-time procedure:
P1:= subs(
    {
        _rp1= 'c', _rp1d= 7, 
        _op2= 'd', _op2t= ':-posint', _op2d= 3, 
        _kp3= 'e', _kp3d= 5,
        _kp4= 'f', _kp4t= ':-integer', _kp4d= 0.5
    },
    eval(P)
);

P1:= proc(c:= 7, d::posint:= 3, {e:= 5, f::integer:= .5})
    c, d, e, f 
end proc

P1();
                          7, 3, 5, 0.5

P1(':-e'= 3, ':-f'= 7);
                           7, 3, 3, 7

Notes:

  • The eval is necessary because of the last name evaluation property of procedures (see ?last_name_eval).
  • As always, keywords must be global symbols, but any issues related to competing local variables or assigned values can be resolved with ':-...(as is usual).
  • Substitutions are effective in both the header and body of the procedure; those in the body can involve names local to the procedure.
  • Just like ordinary subs, there can be multiple sets or lists of substitutions, which are done sequentially (see ?subs).
  • You can substitute procedures for the symbols. In particular, if you already have a partial solution to your problem constructed with codegen, you can use it.
  • You can use subs[eval] instead of subs. This can be particularly effective when combined with the point immediately above.

Method 2: The FromInert method: This method can be unwieldy due to every substructure needing to be expanded to an inert form down to its minutest detail, but on the other hand, it does handle every possible case. FromInert allows for the dynamic creation of any Maple structure including procedures and modules. You need to learn the "language" that FromInert uses, which is where its inverse command ToInert helps. For any structure (other than a "naked" expression sequence),

S1:= FromInert(ToInert(eval(S)))

will return a copy of S. So, to use this for dynamic procedure creation, you start with a template and expand it with ToInert. Here, I'm using the same template as with the subs method:

InertP:= ToInert(eval(P));
InertP:= _Inert_PROC(_Inert_PARAMSEQ(
  _Inert_ASSIGN(_Inert_NAME("_rp1"), _Inert_NAME("_rp1d")), 
  _Inert_ASSIGN(
  _Inert_DCOLON(_Inert_NAME("_op2"), _Inert_NAME("_op2t")), 
  _Inert_NAME("_op2d")), _Inert_SET(_Inert_EXPSEQ(
  _Inert_ASSIGN(_Inert_NAME("_kp3"), _Inert_NAME("_kp3d")), 
  _Inert_ASSIGN(
  _Inert_DCOLON(_Inert_NAME("_kp4"), _Inert_NAME("_kp4t")), 
  _Inert_NAME("_kp4d"))))), _Inert_LOCALSEQ(), 
  _Inert_OPTIONSEQ(), _Inert_EXPSEQ(), _Inert_STATSEQ(
  _Inert_EXPSEQ(_Inert_PARAM(1), _Inert_PARAM(2), 
  _Inert_PARAM(3), _Inert_PARAM(4))), _Inert_DESCRIPTIONSEQ(), 
  _Inert_GLOBALSEQ(), _Inert_LEXICALSEQ(), 
  _Inert_EOP(_Inert_EXPSEQ()))

So, you manipulate this structure using ordinary Maple commands (such as subsindets, and subsindets) and then reassemble the modified structure with FromInert. Note that this allows you to change the number of parameters, which you can't do by the subs method.

Here's a procedure for it (warning: this procedure is fairly deep):

ShowSteps:= (q::`=`(algebraic), v::{name, function})->
    subsindets(
        (InertForm:-MakeInert@isolate)(
            (value@subsindets[2])(q, numeric, `_%a`, nprintf), v 
        ),
        suffixed(_, numeric),
        parse@substring, 2..-1
    )
:
#Your example. Note that all arithmetic operations (other than =) and function names
#must be entered with %.
q1:= 0.379948 = 2 %/ (1 %+ %exp(-0.004 %* cp)) %- 1;

ShowSteps(q1, cp);
#You'll need to view its output in Maple to appreciate it.

value(%);
                        cp = 199.9994377

 

You only need three changes to make your code run without error:

  1. Change eval to eval[recurse]. The recurse is needed because some of the evaluations use variables on their right sides that appear on the left sides of other evaluations in the same set.
  2. Change (x,y,t)-> test to unapply(test, [x,y,t]). The arrow -> cannot be used to make variables that appear indirectly on its right side match the parameters on its left side.
  3. Change plot to plot3d. The former command is strictly for 2D plots.

Use the try command to trap the "time expired" error and break the loop, like this:

for k do
    try timelimit(10, MyProc(k))
    catch "time expired": break
    end try
od:

 

First 73 74 75 76 77 78 79 Last Page 75 of 395