MaplePrimes Posts

MaplePrimes Posts are for sharing your experiences, techniques and opinions about Maple, MapleSim and related products, as well as general interests in math and computing.

Latest Post
  • Latest Posts Feed
  • Philip Yasskin, a long-time Maple user and professor at Texas A&M University is passionate about getting young people engaged in mathematics. One of his programs is SEE-Math: a two-week summer day camp for gifted middle school children interested in math. Maplesoft has been a long-standing supporter of SEE-Math, providing software and prizes for the campers.

    A major project in SEE-Math is developing computer animations using Maple. Students spend their time creating various animations, in hopes of taking the top prize at the end of the workshop. A slew of animations are submitted, some with pop-culture references, elaborate plot lines, and incredible detail. The top animations take home prizes, while all animations from that year are featured on the SEE-Math website.

    Maplesoft proudly sponsors this event, and many like it, to promote interest in STEM education. To see all of the animations from this year’s SEE-Math camp, please visit: http://see-math.math.tamu.edu/2015/. You can find the animations listed under “Euler,” “Godel,” “Noether,” and “Ramanujan,” found halfway down the page.

    On this week I asked Maplesoft Customer Service for help. Here is our correspondence
    (Only the purchase code and e-mail addresses are censored. PS. Also the last name of Kari was deleted by Bryon Thur on 28.08.2015.).
    I think this is of interest for many Maple users. I have got some experience contacting
    with Kaspersky Antivirus (They helped me by the use of indices of my comp.) and ABBYYLingvo
    (They helped to install an ABBYYLingvo vocabulary on my phone.) so I can compare and
    make conclusions.


    From:
    Sent: August-15-15 4:44 AM
    To: Maplesoft Customer Service
    Subject: Customer Service Request: (Web) Installation questions

    Hello,
    After upgrading my Windows 7 HB 32-bit to Windows 10 I cannot uninstall my Maple 16 PE.
     It cannot be uninstalled by neither Start/Parameters/System/Applications nor
    Uninstall in C/ProgramFiles/Maple 16.
    The Uninstall option is not seen in Maple 16 as application.
    Also the overinstallation of Maple 16 does not work.
    Waiting for your feedback.
    Sincerely,
    Markiyan Hirnyk
    ------------------------------------------------------------------------------------------------------

    Dear Markiyan Hirnyk,

    Thank you for contacting Maplesoft.

    Maple 16 is not officially supported on Windows 10 but I have added an activation to your
    existing Maple 16 Personal Edition purchase code: XXXXXXXXXXXXXXXX to see if reactivating
    your license fixes the issue.  If reactivating doesn't give you access to Maple 16, please
    send me the exact wording of any error messages that you receive so that I can send
     the information to our Technical Support Team so that they can investigate further.

    Kind regards,

    Kari
    Maplesoft
    Customer Service
    -----------------------------------------------------------------------------------------------------------
    Hello Kari,
    Unfortunately, neither the  reactivation of my Maple 16 PE by XXXXXXXXXX nor its uninstallation
    do not succeed for me. See the error communications in the attached screens (both in one file) screens_1_2.docx.
    It should be noticed that Maple V Release 4 works on Windows 10 of my comp without any problems.
    Regards,
    Markiyan Hirnyk

    --------------------------------------------------------------------------------------------------------
    Hi Markiyan Hirnyk,

    Thanks for your response.
    I am forwarding your information to our Technical Support Team.
    A representative will contact you soon.
    Kind regards, Kari
    Maplesoft Customer Service
    ------------------------------------------------------------------------------------------------------------
    Hello Markiyan,

    This error is usually caused by a Windows permissions setting. To fix this, please do the following:

    1. Ensure that all Maple programs are completely closed.
    2. Click on your Start Menu and go to the 'Programs' > 'Maple 16' > 'Tools' folder.
    3. Right click the 'Activate Maple' icon and choose 'Run as administrator'.
    4. Activate Maple using your purchase code and this should fix your issue.

    Please let me know if you continue to experience any troubles.

    Regards,

    Chris
    Technical Support Analyst
    --------------------------------------------------------------------------------------------------------
    Hello Chris,
    Following your directions, I have just reactivated Maple 16, but my problem is not solved.
    To shed light on the situation, my Maple 16.02 works properly,
    but I cannot uninstall it after upgrading to Windows 10 Home 32bit.
    See the screens in the attached file screen.docx .
    Regards,
    Markiyan Hirnyk
    ----------------------------------------------------------------------------------------
    Hello Markiyan,
    If you are seeing error messages about Maple still being open,
     I would suggest you try to restart your PC and then attempt the uninstall again to ensure
    that you do not have any lingering Maple programs running. Please let me know
    if you still see this message after restarting.

    Regards,

    Chris
    Technical Support Analyst
    ------------------------------------------------------------------------------------------------------
    Hello Chris,
    This does not help too. My guess is execution failure when Maple 16 was installing.
    Because of that reason the Maple 16 installer did not create Maple uninstaller in my Maple 16.
     See the attached screen of the uninstall folder in C:/ Program Files/ screen_3.docx.
    Regards,
    Markiyan Hirnyk
    --------------------------------------------------------------------------------------------------------
    Hello,
    If you think that you have a corrupted installation, I recommend that you reinstall Maple
    using the new Maple 16.02 installer link provided below.
     This version of the installer was created to get around the Windows 8 installation issues and
     may be of help to you in Windows 10 as well, though again please be aware that
     we do not officially support Windows 10 yet.
    Here are the steps to reinstall Maple:
    1. Click on the Start Menu > Control Panel > Programs and Features ( or Add/Remove Programs).
     Find ‘Maple 16’ in the list and uninstall it. If this is not possible, move on to the next step and continue.
    2. Restart your computer.
    3. Click on the Start Menu > Computer > Local Disk C: > Program Files.
    If there is a folder here called ‘Maple 16’, please delete it.
    4. Download the installer for Maple 16 from the following link:
            http://www.maplesoft.com/downloads/?d=C75DEBEC838C08BB1DCCED0440B49503&pr=Maple16
    5. Make sure to download the correct version for your operating system, i.e. Windows version and 32 or 64-bit.
    6. Install Maple by right clicking the installation file and choosing ‘Run as administrator’.
    I hope that this helps to resolve the issues that you’re having and if it does not,
    contact us and we can further investigate for you.

    Regards,
    Chris
    Technical Support Analyst
    --------------------------------------------------------------------------------------------------------
     Hi Chris,
    My problem with Maple 16 is solved. I completely uninstalled it by Uninstall Tool 3.4, not using brute force. After that I installed Maple 16 by the distributive suggested by you. That's all right.
    Regards,
    Markiyan Hirnyk
    ---------------------------------------------------------------------------------------------------------
    Alright, that is good to hear. Please let us know if you run into any further issues with your installation.
    Regards,
    Chris
    Technical Support Analyst

    In addition to providing access to powerful tools for mathematical computation, Maple has been designed to help you work quickly and efficiently. Here are 10 useful short-cuts when working with Maple:

    1. Use F5 to switch between Text and 2D Math input modes in Maple.

    2. Use F2 (Control+? for Macintosh) to quickly bring up Maple Help information for anything that you have typed in your document.

    3.  Automatic Command Completion can be used when you don't want to type in the full name of a Maple command. To use, begin typing the first few letters of the command name, and press CTRL+Space (Esc or Command+Shift+Space for Macintosh, CTRL+Shift+Space for Linux).  A list of possible completions will display; click the one you want.

     

    4. The Shift+Enter key combination lets you continue entering math or commands on a new line without executing that line. 

    5. If you want more than a single command to be executed at once, you must separate them with a semi-colon or colon.

    6. When you click inside a set of commands in Math mode, the dash line indicates the boundaries of the input region; all commands in this region will execute together in sequence.

    7. To increase the size of a piecewise function, add a new row.  Place the cursor on the last row, and press CTRL+Shift+R (Command+Shift+R for Macintosh). These shortcut keys also work to add rows to matrices.

    8. An easy way to insert a Greek letter is to first press CTRL+Shift+G (Command+Shift+G for the Macintosh). The next letter typed will appear in Greek.

    9. Sometimes you may want to insert symbols above or below another character, for example, to enter a vector arrow. To insert a symbol above (called "overscript"), press CTRL+Shift+["] (Command+Shift+["] for Macintosh) and then type in your symbol (or insert it from a palette).

    For example, typing "x" then holding down CTRL+Shift and pressing ["] allows you to insert a symbol above the x, such as 

    10. Compute or recompute the entire Maple worksheet when you have changed expressions that affect subsequent Maple commands.  Press Ctrl + Shift + Enter (Command + Shift + Enter in Macintosh) or click the execute worksheet icon. 

    Are there any short-cuts that you would add to this list?

    In this work the theme of vector analysis shown from a computational point of view; this being a very important role in the engineering component; in civil and mechanical special it is why, using the scientific software Maple develops interactive solutions for long processes through MapleCloud calculations. At present the majority of professors / researchers perform static classes open source leaves; so that our students learn and memorize commands, thus generating more time learning in the area. Loading Bookseller VectorCalculus develop topics: vector algebra, differential operators, conservative fields, etc. Maplesoft making processes provide immediate calculations long operation Embedded Components displayed in line with MapleNet integrations. Today our future engineers to design solutions and will be launched in the cloud thus being a process with global qualification in the specialty. Significantly Maple is a scientific software which allows the researcher to design their own innovations and not use themes for their manufacturers.

     

    III_CRF_2015.pdf

    CRF_2015.mw

     

    L.AraujoC.

     

     

    knight's tour is a sequence of moves of a knight on a chessboard such that the knight visits every square only once. This problem mentioned in  page of tasks still without  Maple implementation. 

    The post presents the implementation of this task in Maple. Required parameter of the procedure (named  KnightTour)  is  address  - the address of the initial square in the algebraic notation. The second parameter  opt  is optional parameter:

    1) if  opt is sequence  (by default) then  the procedure returns the sequence of moves of the knight in the usual algebraic notation,

    2)  if  opt is diagram   then  the procedure returns the plot of moves of the knight and  sequentially numbers all the visited squares,

    3) if  opt is animation  then  the procedure returns an animation of moves of the knight.

    In the procedure is used a solution with maximum symmetry by George Jelliss, http://www.mayhematics.com/t/8f.htm

     

    Code of the procedure:

    KnightTour := proc(address::symbol, opt::symbol := sequence)

    local L, n, L1, k, i, j, squares, border, chessboard, letters, digits, L2, L3, Tour, T, F;

    uses ListTools, plottools, plots;

    L := [a1, b3, d2, c4, a5, b7, d8, e6, d4, b5, c7, a8, b6, c8, a7, c6, b8, a6, b4, d5, e3, d1, b2, a4, c5, d7, f8, h7, f6, g8, h6, f7, h8, g6, e7, f5, h4, g2, e1, d3, e5, g4, f2, h1, g3, f1, h2, f3, g1, h3, g5, e4, d6, e8, g7, h5, f4, e2, c1, a2, c3, b1, a3, c2];

    n := Search(address, L);

    L1 := [L[n .. 64][], L[1 .. n-1][]];

    if opt = sequence then return L1[] fi;

    k := 0;

    for i to 8 do

    for j from `if`(type(i, odd), 1, 2) by 2 to 8 do

    k := k+1;

    squares[k] := polygon([[i-1/2, j-1/2], [i-1/2, j+1/2], [i+1/2, j+1/2], [i+1/2, j-1/2]], color = grey);

    od;  od;

    squares := convert(squares, list);

    border := curve([[1/2, 1/2], [1/2, 17/2], [17/2, 17/2], [17/2, 1/2], [1/2, 1/2]], color = black, thickness = 4);

    chessboard := display(squares, border);

    letters := [a, b, c, d, e, f, g, h];

    digits := [$ 1 .. 8];

    L2 := convert~(L1, string);

    L3 := subs(letters=~digits, map(t->[parse(t[1]), parse(t[2])], L2));

    Tour := curve(L3, color = red, thickness = 3);

    T := textplot([seq([op(L3[i]), i], i = 1 .. 64)], align = above, font = [times, bold, 16]);

    if opt = diagram then return display(chessboard, Tour, T, axes = none) fi;

    F := seq(display(chessboard, curve(L3[1 .. s], color = red, thickness = 3), textplot([seq([op(L3[i]), i], i = 1 .. s)], align = above, font = [times, bold, 16])), s = 1 .. 64);

    display(seq(F[i]$5, i = 1 .. 64), insequence = true, axes = none);

    end proc:

     

     Examples of use:

    KnightTour(f3);

    KnightTour(f3, diagram);

     

     

    KnightTour(f3, animation);

                                     

     

     

     KnightTour.mw

    As an educator, you surely know that giving students more involved problems in an online assessment tool provides challenges, both for students and instructors. Only marking the final answer doesn’t necessarily provide an understanding of the student’s capabilities, and it penalizes students that make a small mistake in part of their solution.

    This webinar will demonstrate how to create separate questions with multiple steps that can be linked or chained together. Question chaining allows instructors to mark subsequent questions based on the correct answer or the answer provided by students in previous parts.

    To join us for the live presentation, please click here to register.

    Dear friends,

    some time ago I shared a story here on the use of Maple to compute the cycle index of the induced action on the edges of an ordinary graph of the symmetric group permuting the vertices and the use of the Polya Enumeration Theorem to count non-isomorphic graphs by the number of edges. It can be found at the following Mapleprimes link.

    I am writing today to alert you to another simple Maple program that is closely related and demonstrates Maple's capability to implement concepts from group theory and Polya enumeration. This link at Math.Stackexchange.com shows how to use the cycle index of the induced action by the symmetric group permuting vertices on the edges of a multigraph that includes loops to count set partitions of multisets containing two instances of n distinct types of items. The sequence that corresponds to these set partitions is OEIS A020555 where it is pointed out that we can equivalently count multigraphs with n labeled i.e. distinct edges where the vertices of the graph represent the multisets of the multiset partition and are connected by an edge k if the two instances of the value k are included in the sets represented by the two vertices that constitute the edge. The problem then reduces to a simple substitution into the aforementioned cycle index of a polynomial representing the set of labels on an edge including no labels on an edge that is not included.

    This computation presents a remarkable simplicity while also implementing a non-trivial application of Polya counting. It is hoped that MaplePrimes users will enjoy reading this program, possibly profit from some of the techniques employed and be motivated to use Maple in their work on combinatorics problems.

    Best regards,

    Marko Riedel

    In this paper of presents under a totally modern sound environment dynamics; using embedded components that gives us the Cybernet Company through its product Maple 2015. Using classical techniques vector equations describe the particle, particle system and solid bodies. We note that the solutions o ered by this software motivate students civil and mechanical engineering to nd optimal answers. Integrating algorithms own programming language and solid mechanics using buttons we relate the movement of translation and rotation with reference to its center of mass.
    Choosing envelopes graphical methods, functional programming and mathematical computer display modeling reached alternatives to achieve the next generation of engineers. Therefore this work show that the use of
    embedded components allow us to merge the traditional and the computer; It means that all these equations using physical and propose viable criteria we perform in a dynamic sheet; which they have a number of components; then generate simulations with real objects.

    Congreso COMAP 2015.pdf

    Study of the Dynamics of the Solid with Embedded Components in Civil Engineering with Maplesoft.mw

    (in spanish)

    L.AraujoC.

     

    Symbolic sequences enter in various formulations in mathematics. This post is about a related new subpackage, Sequences, within the MathematicalFunctions package, available for download in Maplesoft's R&D page for Mathematical Functions and Differential Equations (currently bundled with updates to the Physics package).

     

    Perhaps the most typical cases of symbolic sequences are:

     

    1) A sequence of numbers - say from n to m - frequently displayed as

    n, `...`, m

     

    2) A sequence of one object, say a, repeated say p times, frequently displayed as

     "((a,`...`,a))"

    3) A more general sequence, as in 1), but of different objects and not necessarily numbers, frequently displayed as

    a[n], `...`, a[m]

    or likewise a sequence of functions

    f(n), `...`, f(m)

    In all these cases, of course, none of n, m, or p are known: they are just symbols, or algebraic expressions, representing integer values.

     

    These most typical cases of symbolic sequences have been implemented in Maple since day 1 using the `$` operator. Cases 1), 2) and 3) above are respectively entered as `$`(n .. m), `$`(a, p), and `$`(a[i], i = n .. m) or "`$`(f(i), i = n .. m)." To have computer algebra representations for all these symbolic sequences is something wonderful, I would say unique in Maple.

    Until recently, however, the typesetting of these symbolic sequences was frankly poor, input like `$`(a[i], i = n .. m) or ``$\``(a, p) just being echoed in the display. More relevant: too little could be done with these objects; the rest of Maple didn't know how to add, multiply, differentiate or map an operation over the elements of the sequence, nor for instance count the sequence's number of elements.

     

    All this has now been implemented.  What follows is a brief illustration.

    restart

    First of all, now these three types of sequences have textbook-like typesetting:

    `$`(n .. m)

    `$`(n .. m)

    (1)

    `$`(a, p)

    `$`(a, p)

    (2)

    For the above, a$p works the same way

    `$`(a[i], i = n .. m)

    `$`(a[i], i = n .. m)

    (3)

    Moreover, this now permits textbook display of mathematical functions that depend on sequences of paramateters, for example:

    hypergeom([`$`(a[i], i = 1 .. p)], [`$`(b[i], i = 1 .. q)], z)

    hypergeom([`$`(a[i], i = 1 .. p)], [`$`(b[i], i = 1 .. q)], z)

    (4)

    IncompleteBellB(n, k, `$`(factorial(j), j = 1 .. n-k+1))

    IncompleteBellB(n, k, `$`(factorial(j), j = 1 .. n-k+1))

    (5)

    More interestingly, these new developments now permit differentiating these functions even when their arguments are symbolic sequences, and displaying the result as in textbooks, with copy and paste working properly, for instance

    (%diff = diff)(hypergeom([`$`(a[i], i = 1 .. p)], [`$`(b[i], i = 1 .. q)], z), z)

    %diff(hypergeom([`$`(a[i], i = 1 .. p)], [`$`(b[i], i = 1 .. q)], z), z) = (product(a[i], i = 1 .. p))*hypergeom([`$`(a[i]+1, i = 1 .. p)], [`$`(b[i]+1, i = 1 .. q)], z)/(product(b[i], i = 1 .. q))

    (6)

    It is very interesting how much this enhances the representation capabilities; to mention but one, this makes 100% possible the implementation of the Faa-di-Bruno  formula for the nth symbolic derivative of composite functions (more on this in a post to follow this one).

    But the bread-and-butter first: the new package for handling sequences is

    with(MathematicalFunctions:-Sequences)

    [Add, Differentiate, Map, Multiply, Nops]

    (7)

    The five commands that got loaded do what their name tells. Consider for instance the first kind of sequences mentione above, i.e

    `$`(n .. m)

    `$`(n .. m)

    (8)

    Check what is behind this nice typesetting

    lprint(`$`(n .. m))

    `$`(n .. m)

     

    All OK. How many operands (an abstract version of Maple's nops  command):

    Nops(`$`(n .. m))

    m-n+1

    (9)

    That was easy, ok. Add the sequence

    Add(`$`(n .. m))

    (1/2)*(m-n+1)*(n+m)

    (10)

    Multiply the sequence

    Multiply(`$`(n .. m))

    factorial(m)/factorial(n-1)

    (11)

    Map an operation over the elements of the sequence

    Map(f, `$`(n .. m))

    `$`(f(j), j = n .. m)

    (12)

    lprint(`$`(f(j), j = n .. m))

    `$`(f(j), j = n .. m)

     

    Map works as map, i.e. you can map extra arguments as well

    MathematicalFunctions:-Sequences:-Map(Int, `$`(n .. m), x)

    `$`(Int(j, x), j = n .. m)

    (13)

    All this works the same way with symbolic sequences of forms "((a,`...`,a))" , and a[n], `...`, a[m]. For example:

    `$`(a, p)

    `$`(a, p)

    (14)

    lprint(`$`(a, p))

    `$`(a, p)

     

    MathematicalFunctions:-Sequences:-Nops(`$`(a, p))

    p

    (15)

    Add(`$`(a, p))

    a*p

    (16)

    Multiply(`$`(a, p))

    a^p

    (17)

    Differentation also works

    Differentiate(`$`(a, p), a)

    `$`(1, p)

    (18)

    MathematicalFunctions:-Sequences:-Map(f, `$`(a, p))

    `$`(f(a), p)

    (19)

    MathematicalFunctions:-Sequences:-Differentiate(`$`(f(a), p), a)

    `$`(diff(f(a), a), p)

    (20)

    For a symbolic sequence of type 3)

    `$`(a[i], i = n .. m)

    `$`(a[i], i = n .. m)

    (21)

    MathematicalFunctions:-Sequences:-Nops(`$`(a[i], i = n .. m))

    m-n+1

    (22)

    Add(`$`(a[i], i = n .. m))

    sum(a[i], i = n .. m)

    (23)

    Multiply(`$`(a[i], i = n .. m))

    product(a[i], i = n .. m)

    (24)

    The following is nontrivial: differentiating the sequence a[n], `...`, a[m], with respect to a[k] should return 1 when n = k (i.e the running index has the value k), and 0 otherwise, and the same regarding m and k. That is how it works now:

    Differentiate(`$`(a[i], i = n .. m), a[k])

    `$`(piecewise(k = i, 1, 0), i = n .. m)

    (25)

    lprint(`$`(piecewise(k = i, 1, 0), i = n .. m))

    `$`(piecewise(k = i, 1, 0), i = n .. m)

     

    MathematicalFunctions:-Sequences:-Map(f, `$`(a[i], i = n .. m))

    `$`(f(a[i]), i = n .. m)

    (26)

    Differentiate(`$`(f(a[i]), i = n .. m), a[k])

    `$`((diff(f(a[i]), a[i]))*piecewise(k = i, 1, 0), i = n .. m)

    (27)

    lprint(`$`((diff(f(a[i]), a[i]))*piecewise(k = i, 1, 0), i = n .. m))

    `$`((diff(f(a[i]), a[i]))*piecewise(k = i, 1, 0), i = n .. m)

     

     

    And that is it. Summarizing: in addition to the former implementation of symbolic sequences, we now have textbook-like typesetting for them, and more important: Add, Multiply, Differentiate, Map and Nops. :)

     

    The first large application we have been working on taking advantage of this is symbolic differentiation, with very nice results; I will see to summarize them in a post to follow in a couple of days.

     

    Download MathematicalFunctionsSequences.mw

    Edgardo S. Cheb-Terrab
    Physics, Differential Equations and Mathematical Functions, Maplesoft

    For the past thirty years, I have used several mathematical packages for problem solving and graphing. It all started with spreadsheet software that really helped speedup calculations compared to calculators. As many people do, once I had one tool I then started looking for another that would offer even more capabilities and features. I tested several of the very early math software but none really did all that I wanted until I came across Maple while I was working at SPAR Aerospace in Canada. For me, the rest is history. As long as I had a copy of Maple, it was all that I needed.

    On occasions when I did not have a copy of this amazing software, I resorted to spreadsheets once more to complete fairly large and complex projects involving large databases and large numbers of calculations, especially when performing What-If scenarios. One distinct disadvantage of using a spreadsheet was the cryptic form of equation writing. I had to divide one long equation into several sections in different cells and then add them all up, which clearly is not good for documentation of the calculations. It is also very confusing for other engineers to know what that equation is or what it does. The development of the full engineering spreadsheet took months to complete, debug and verify. During this process, when I had errors, it was often very difficult to track exactly where the problem was, making the debugging process time consuming and sometimes very frustrating.

    Having worked with Maple before, I remembered how easy it was to enter equations in a very familiar, readable math format. The real power of this software is that it allows you to write the equation(s) anyway you like and solve for any given parameter, unlike spreadsheets where you have to solve the problem first, by hand, for the parameter you want and then get the spreadsheet to calculate the value. I remember one time a few years ago when I wrote nine or ten simultaneous differential equations all in symbolic form and asked Maple to calculate certain parameters in a fully symbolic form. To my utmost disbelief, the answer came back within few minutes. With results in hand, I was able to quickly finish my research, and the results were published at PCIM Europe 2005 in “Distributed Gate ESR and its Effect on Shoot Through Performance at the Die Level”. I would never have gotten the results I needed if I was using a spreadsheet.

    Even with much simpler systems of equations, finding solutions with a paper and pencil was never an easy task for me. It took a very long time, and even then there was no guarantee that I did not make copying errors, accidentally leave out a term, or make a calculation error. After I found the correct solution, I then had the problem of plotting the results, which I often needed in 3-D. Plotting allowed much deeper insights into the interdependency of all the parameters and made it easy for me to concentrate on the important ones without wasting any time. I was very happy when I could pass all these tasks onto Maple, which could do them much faster and more reliably then I ever could. Maple is a software that allows me to go beyond routine engineering calculations and gives me the tools to reach levels of insight and understanding that were completely out of reach of the average engineer until a few years ago.

    For the record, I have no business affiliations with Maplesoft. I’m writing this article because Maple makes such a difference in my work that I feel it is important to share my experiences so other engineers can get the same benefits.

           Calculation of RSCR mechanism as a  solution to underdetermined system of nonlinear equations.  

     

    RSCR.mw 

    https://vk.com/doc242471809_376439263
    https://vk.com/doc242471809_408704758

    RCCC mechanism
    https://vk.com/doc242471809_375452868

    Here we develop the factoring in common factor, simple and complete square blade, plus simple equation systems with graphic and design, and graphic solution of the quadratic equation using components in maple 2015.

     

    Factorizacion.mw

    (in spanish)

    L.AraujoC.

     

     

     

    Solving DAEs in Maple

    As I had mentioned in many posts, Maple cannot solve nonlinear DAEs. As of today (Maple 2015), given a system of index 1 DAE dy/dt = f(y,z); 0 = g(y,z), Maple “extends” g(y,z) to get dz/dt = g1(y,z). So, any given index 1 DAE is converted to a system of ODEs dy/dt = f, dz/dt = g1 with the constraint g = 0, before it solves. This is true for all the solvers in Maple despite the wrong claims in the help files. It is also true for MEBDFI, the FORTRAN implementation of which actually solves index 2 DAEs directly. In addition, the initial condition for the algebraic variable has to be consistent. The problem with using fsolve is that one cannot specify tolerance. Often times one has to solve DAEs at lower tolerances (1e-3 or 1e-4), not 1e-15. In addition, one cannot use compile =true for high efficiency.

    The approach in Maple fails for many DAEs of practical interest. In addition, this can give wrong answers. For example, dy/dt = z, y^2+z^2-1=0, with y(0)=0 and z(0)=1 has a meaningful solution only from 0 to Pi/2, Maple’s dsolve/numeric will convert this to dy/dt = z and dz/dt = -y and integrate in time for ever. This is wrong as at t = Pi/2, the system becomes index 2 DAE and there is more than 1 acceptable solution.

    We just recently got our paper accepted that helps Maple's dsolve and many ODE solvers in other languages handle DAEs directly. The approach is rather simple, the index 1 DAE is perturbed as dy/dt = f.S ; -mu.diff(g,t) = g. The switch function is a tanh function which is zero till t = tinit (initialization time). Mu is chosen to be a small number. In addition, lhs of the perturbed equation is simplified using approximate initial guesses as Maple cannot handle non-constant Mass matrix. The paper linked below gives below more details.

    http://depts.washington.edu/maple/pubs/U_Apprach_FULL_DRAFT.pdf  

    Next, I discuss different examples (code attached).

    Example 1: Simple DAE (one ODE and one AE), the proposed approach helps solve this system efficiently without knowing the exact initial condition.

    Hint: The code is given with a semicolon at the end of the most of the statements for educational purposes for new users. Using semicolon after dsolve numeric in classic worksheet crashes the code (Maplesoft folks couldn’t fix this despite my request).

    Example 2:

    This is a nickel battery model (one ODE and one AE). This fails with many solvers unless exact IC is given. The proposed approach works well. In particular, stiff=true, implicit=true option is found to be very efficient. The code given in example 1 can be used to solve example 2 or any other DAEs by just entering ODEs, AEs, ICs in proper places.

    Example 3:

    This is a nonlinear implicit ODE (posted in Mapleprimes earlier by joha, (http://www.mapleprimes.com/questions/203096-Solving-Nonlinear-ODE#answer211682 ). This can be converted to index 1 DAE and solved using the proposed approach.

    Example 4:

    This example was posted earlier by me in Mapleprimes (http://www.mapleprimes.com/posts/149877-ODEs-PDEs-And-Special-Functions) . Don’t try to solve this in Maple using its dsolve numeric solver for N greater than 5 directly. The proposed approach can handle this well. Tuning the perturbation parameters and using compile =true will help in solving this system more efficiently.

    Finally example 1 is solved for different perturbation parameters to show how one can arrive at convergence. Perhaps more sophisticated users and Maplesoft folks can enable this as automatically tuned parameters in dsolve/numeric.

    Note:

    This post should not be viewed as just being critical on Maple. I have wasted a lot of time assuming that Maple does what it claims in its help file. This post is to bring awareness about the difficulty in dealing with DAEs. It is perfectly fine to not have a DAE solver, but when literature or standard methods are cited/claimed, they should be implemented in proper form. I will forever remain a loyal Maple user as it has enabled me to learn different topics efficiently. Note that Maplesim can solve DAEs, but it is a pain to create a Maplesim model/project just for solving a DAE. Also, events is a pain with Maplesim. The reference to Mapleprimes links are missing in the article, but will be added before the final printing stage. The ability of Maple to find analytical Jacobian helps in developing many robust ODE and DAE solvers and I hope to post my own approaches that will solve more complicated systems.

    I strongly encourage testing of the proposed approach and implementation for various educational/research purposes by various users. The proposed approach enables solving lithium-ion and other battery/electrochemical storage models accurately in a robust manner. A disclosure has been filed with the University of Washington to apply for a provisional patent for battery models and Battery Management System for transportation, storage and other applications because of the current commercial interest in this topic (for batteries). In particular, use of this single step avoids intialization issues/(no need to initialize separately) for parameter estimation, state estimation or optimal control of battery models.

     

    Appendix A

    Maple code for Examples 1-4 from "Extending Explicit and Linealry Implicit ODE Solvers for Index-1 DAEs", M. T. Lawder,

    V. Ramadesigan, B. Suthar and V. R. Subramanian, Computers and Chemical Engineering, in press (2015).

    Use y1, y2, etc. for all differential variables and z1, z2, etc. for all algebraic variables

     

    Example 1

    restart;

    with(plots):

    Enter all ODEs in eqode

    eqode:=[diff(y1(t),t)=-y1(t)^2+z1(t)];

    eqode := [diff(y1(t), t) = -y1(t)^2+z1(t)]

    (1)

    Enter all AEs in eqae

    eqae:=[cos(y1(t))-z1(t)^0.5=0];

    eqae := [cos(y1(t))-z1(t)^.5 = 0]

    (2)

    Enter all initial conditions for differential variables in icodes

    icodes:=[y1(0)=0.25];

    icodes := [y1(0) = .25]

    (3)

    Enter all intial conditions for algebraic variables in icaes

    icaes:=[z1(0)=0.8];

    icaes := [z1(0) = .8]

    (4)

    Enter parameters for perturbation value (mu), switch function (q and tint), and runtime (tf)

    pars:=[mu=0.1,q=1000,tint=1,tf=5];

    pars := [mu = .1, q = 1000, tint = 1, tf = 5]

    (5)

    Choose solving method (1 for explicit, 2 for implicit)

    Xexplicit:=2;

    Xexplicit := 2

    (6)

    Standard solver requires IC z(0)=0.938791 or else it will fail

    solx:=dsolve({eqode[1],eqae[1],icodes[1],icaes[1]},numeric):

    Error, (in dsolve/numeric/DAE/checkconstraints) the initial conditions do not satisfy the algebraic constraints
      error = .745e-1, tolerance = .559e-6, constraint = cos(y1(t))-z1(t)^.5000000000000000000000

     

    ff:=subs(pars,1/2+1/2*tanh(q*(t-tint)));

    ff := 1/2+(1/2)*tanh(1000*t-1000)

    (7)

    NODE:=nops(eqode);NAE:=nops(eqae);

    NODE := 1

    NAE := 1

    (8)

    for XX from 1 to NODE do
    EQODE||XX:=lhs(eqode[XX])=rhs(eqode[XX])*ff;
    end do;

    EQODE1 := diff(y1(t), t) = (-y1(t)^2+z1(t))*(1/2+(1/2)*tanh(1000*t-1000))

    (9)

    for XX from 1 to NAE do
    EQAE||XX:=subs(pars,-mu*(diff(rhs(eqae[XX])-lhs(eqae[XX]),t))=rhs(eqae[XX])-lhs(eqae[XX]));
    end do;

    EQAE1 := -.1*sin(y1(t))*(diff(y1(t), t))-0.5e-1*(diff(z1(t), t))/z1(t)^.5 = -cos(y1(t))+z1(t)^.5

    (10)

     

    Dvars1:={seq(diff(z||x(t),t)=D||x,x=1..NAE)};

    Dvars1 := {diff(z1(t), t) = D1}

    (11)

    Dvars2:={seq(rhs(Dvars1[x])=lhs(Dvars1[x]),x=1..NAE)};

    Dvars2 := {D1 = diff(z1(t), t)}

    (12)

    icsn:=seq(subs(y||x(0)=y||x(t),icodes[x]),x=1..NODE),seq(subs(z||x(0)=z||x(t),icaes[x]),x=1..NAE);

    icsn := y1(t) = .25, z1(t) = .8

    (13)

    for j from 1 to NAE do

    EQAEX||j:=subs(Dvars1,eqode,icsn,Dvars2,lhs(EQAE||j))=rhs(EQAE||j);

    end do:

    Sys:={seq(EQODE||x,x=1..NODE),seq(EQAEX||x,x=1..NAE),seq(icodes[x],x=1..NODE),seq(icaes[x],x=1..NAE)};

    Sys := {-0.1824604200e-1-0.5590169945e-1*(diff(z1(t), t)) = -cos(y1(t))+z1(t)^.5, y1(0) = .25, z1(0) = .8, diff(y1(t), t) = (-y1(t)^2+z1(t))*(1/2+(1/2)*tanh(1000*t-1000))}

    (14)

    if Xexplicit=1 then

    sol:=dsolve(Sys,numeric):

    else

    sol:=dsolve(Sys,numeric,stiff=true,implicit=true):
    end if:

     

    for XX from 1 to NODE do
    a||XX:=odeplot(sol,[t,y||XX(t)],0..subs(pars,tf),color=red);
    end do:

    for XX from NODE+1 to NODE+NAE do
    a||XX:=odeplot(sol,[t,z||(XX-NODE)(t)],0..subs(pars,tf),color=blue);
    end do:

    display(seq(a||x,x=1..NODE+NAE),axes=boxed);

     

    End Example 1

     

    Example 2

    restart;

    with(plots):

    eq1:=diff(y1(t),t)=j1*W/F/rho/V;

    eq1 := diff(y1(t), t) = j1*W/(F*rho*V)

    (15)

    eq2:=j1+j2=iapp;

    eq2 := j1+j2 = iapp

    (16)

    j1:=io1*(2*(1-y1(t))*exp((0.5*F/R/T)*(z1(t)-phi1))-2*y1(t)*exp((-0.5*F/R/T)*(z1(t)-phi1)));

    j1 := io1*(2*(1-y1(t))*exp(.5*F*(z1(t)-phi1)/(R*T))-2*y1(t)*exp(-.5*F*(z1(t)-phi1)/(R*T)))

    (17)

    j2:=io2*(exp((F/R/T)*(z1(t)-phi2))-exp((-F/R/T)*(z1(t)-phi2)));

    j2 := io2*(exp(F*(z1(t)-phi2)/(R*T))-exp(-F*(z1(t)-phi2)/(R*T)))

    (18)

    params:={F=96487,R=8.314,T=298.15,phi1=0.420,phi2=0.303,W=92.7,V=1e-5,io1=1e-4,io2=1e-10,iapp=1e-5,rho=3.4};

    params := {F = 96487, R = 8.314, T = 298.15, V = 0.1e-4, W = 92.7, io1 = 0.1e-3, io2 = 0.1e-9, rho = 3.4, iapp = 0.1e-4, phi1 = .420, phi2 = .303}

    (19)

    eqode:=[subs(params,eq1)];

    eqode := [diff(y1(t), t) = 0.5651477584e-2*(1-y1(t))*exp(19.46229155*z1(t)-8.174162450)-0.5651477584e-2*y1(t)*exp(-19.46229155*z1(t)+8.174162450)]

    (20)

    eqae:=[subs(params,eq2)];

    eqae := [0.2e-3*(1-y1(t))*exp(19.46229155*z1(t)-8.174162450)-0.2e-3*y1(t)*exp(-19.46229155*z1(t)+8.174162450)+0.1e-9*exp(38.92458310*z1(t)-11.79414868)-0.1e-9*exp(-38.92458310*z1(t)+11.79414868) = 0.1e-4]

    (21)

    icodes:=[y1(0)=0.05];

    icodes := [y1(0) = 0.5e-1]

    (22)

    icaes:=[z1(0)=0.7];

    icaes := [z1(0) = .7]

    (23)

    solx:=dsolve({eqode[1],eqae[1],icodes[1],icaes[1]},type=numeric):

    Error, (in dsolve/numeric/DAE/checkconstraints) the initial conditions do not satisfy the algebraic constraints
      error = .447e9, tolerance = .880e4, constraint = -2000000*(-1+y1(t))*exp(19.46229155000000000000*z1(t)-8.174162450000000000000)-2000000*y1(t)*exp(-19.46229155000000000000*z1(t)+8.174162450000000000000)+exp(38.92458310000000000000*z1(t)-11.79414868000000000000)-exp(-38.92458310000000000000*z1(t)+11.79414868000000000000)-100000

     

    pars:=[mu=0.00001,q=1000,tint=1,tf=5001];

    pars := [mu = 0.1e-4, q = 1000, tint = 1, tf = 5001]

    (24)

    Xexplicit:=2;

    Xexplicit := 2

    (25)

    ff:=subs(pars,1/2+1/2*tanh(q*(t-tint)));

    ff := 1/2+(1/2)*tanh(1000*t-1000)

    (26)

    NODE:=nops(eqode):NAE:=nops(eqae);

    NAE := 1

    (27)

    for XX from 1 to NODE do
    EQODE||XX:=lhs(eqode[XX])=rhs(eqode[XX])*ff:
    end do;

    EQODE1 := diff(y1(t), t) = (0.5651477584e-2*(1-y1(t))*exp(19.46229155*z1(t)-8.174162450)-0.5651477584e-2*y1(t)*exp(-19.46229155*z1(t)+8.174162450))*(1/2+(1/2)*tanh(1000*t-1000))

    (28)

    for XX from 1 to NAE do
    EQAE||XX:=subs(pars,-mu*(diff(rhs(eqae[XX])-lhs(eqae[XX]),t))=rhs(eqae[XX])-lhs(eqae[XX]));
    end do;

    EQAE1 := -0.2e-8*(diff(y1(t), t))*exp(19.46229155*z1(t)-8.174162450)+0.3892458310e-7*(1-y1(t))*(diff(z1(t), t))*exp(19.46229155*z1(t)-8.174162450)-0.2e-8*(diff(y1(t), t))*exp(-19.46229155*z1(t)+8.174162450)+0.3892458310e-7*y1(t)*(diff(z1(t), t))*exp(-19.46229155*z1(t)+8.174162450)+0.3892458310e-13*(diff(z1(t), t))*exp(38.92458310*z1(t)-11.79414868)+0.3892458310e-13*(diff(z1(t), t))*exp(-38.92458310*z1(t)+11.79414868) = 0.1e-4-0.2e-3*(1-y1(t))*exp(19.46229155*z1(t)-8.174162450)+0.2e-3*y1(t)*exp(-19.46229155*z1(t)+8.174162450)-0.1e-9*exp(38.92458310*z1(t)-11.79414868)+0.1e-9*exp(-38.92458310*z1(t)+11.79414868)

    (29)

    Dvars1:={seq(diff(z||x(t),t)=D||x,x=1..NAE)};

    Dvars1 := {diff(z1(t), t) = D1}

    (30)

    Dvars2:={seq(rhs(Dvars1[x])=lhs(Dvars1[x]),x=1..NAE)};

    Dvars2 := {D1 = diff(z1(t), t)}

    (31)

    icsn:=seq(subs(y||x(0)=y||x(t),icodes[x]),x=1..NODE),seq(subs(z||x(0)=z||x(t),icaes[x]),x=1..NAE);

    icsn := y1(t) = 0.5e-1, z1(t) = .7

    (32)

    for j from 1 to NAE do

    EQAEX||j:=subs(Dvars1,eqode,icsn,Dvars2,lhs(EQAE||j))=rhs(EQAE||j);

    end do;

    EQAEX1 := -0.2e-8*(0.5368903705e-2*exp(5.449441630)-0.2825738792e-3*exp(-5.449441630))*exp(5.449441630)+0.3697835394e-7*(diff(z1(t), t))*exp(5.449441630)-0.2e-8*(0.5368903705e-2*exp(5.449441630)-0.2825738792e-3*exp(-5.449441630))*exp(-5.449441630)+0.1946229155e-8*(diff(z1(t), t))*exp(-5.449441630)+0.3892458310e-13*(diff(z1(t), t))*exp(15.45305949)+0.3892458310e-13*(diff(z1(t), t))*exp(-15.45305949) = 0.1e-4-0.2e-3*(1-y1(t))*exp(19.46229155*z1(t)-8.174162450)+0.2e-3*y1(t)*exp(-19.46229155*z1(t)+8.174162450)-0.1e-9*exp(38.92458310*z1(t)-11.79414868)+0.1e-9*exp(-38.92458310*z1(t)+11.79414868)

    (33)

    Sys:={seq(EQODE||x,x=1..NODE),seq(EQAEX||x,x=1..NAE),seq(icodes[x],x=1..NODE),seq(icaes[x],x=1..NAE)};

    Sys := {-0.5810962488e-6+0.8802389238e-5*(diff(z1(t), t)) = 0.1e-4-0.2e-3*(1-y1(t))*exp(19.46229155*z1(t)-8.174162450)+0.2e-3*y1(t)*exp(-19.46229155*z1(t)+8.174162450)-0.1e-9*exp(38.92458310*z1(t)-11.79414868)+0.1e-9*exp(-38.92458310*z1(t)+11.79414868), y1(0) = 0.5e-1, z1(0) = .7, diff(y1(t), t) = (0.5651477584e-2*(1-y1(t))*exp(19.46229155*z1(t)-8.174162450)-0.5651477584e-2*y1(t)*exp(-19.46229155*z1(t)+8.174162450))*(1/2+(1/2)*tanh(1000*t-1000))}

    (34)

    if Xexplicit=1 then

    sol:=dsolve(Sys,numeric,maxfun=0):

    else

    sol:=dsolve(Sys,numeric,stiff=true,implicit=true,maxfun=0):

    end if:

     

    for XX from 1 to NODE do
    a||XX:=odeplot(sol,[t,y||XX(t)],0..subs(pars,tf),color=red);
    end do:

    for XX from NODE+1 to NODE+NAE do
    a||XX:=odeplot(sol,[t,z||(XX-NODE)(t)],0..subs(pars,tf),color=blue);
    end do:

    b1:=odeplot(sol,[t,y1(t)],0..1,color=red):
    b2:=odeplot(sol,[t,z1(t)],0..1,color=blue):

    display(b1,b2,axes=boxed);

     

    display(seq(a||x,x=1..NODE+NAE),axes=boxed);

     

    End Example 2

     

    Example 3

    restart;

    with(plots):

    eq1:=diff(y1(t),t)^2+diff(y1(t),t)*(y1(t)+1)+y1(t)=cos(diff(y1(t),t));

    eq1 := (diff(y1(t), t))^2+(diff(y1(t), t))*(y1(t)+1)+y1(t) = cos(diff(y1(t), t))

    (35)

    solx:=dsolve({eq1,y1(0)=0},numeric):

    Error, (in dsolve/numeric/make_proc) Could not convert to an explicit first order system due to 'RootOf'

     

    eqode:=[diff(y1(t),t)=z1(t)];

    eqode := [diff(y1(t), t) = z1(t)]

    (36)

    eqae:=[subs(eqode,eq1)];

    eqae := [z1(t)^2+z1(t)*(y1(t)+1)+y1(t) = cos(z1(t))]

    (37)

    icodes:=[y1(0)=0.0];

    icodes := [y1(0) = 0.]

    (38)

    icaes:=[z1(0)=0.0];

    icaes := [z1(0) = 0.]

    (39)

    pars:=[mu=0.1,q=1000,tint=1,tf=4];

    pars := [mu = .1, q = 1000, tint = 1, tf = 4]

    (40)

    Xexplicit:=2;

    Xexplicit := 2

    (41)

    ff:=subs(pars,1/2+1/2*tanh(q*(t-tint)));

    ff := 1/2+(1/2)*tanh(1000*t-1000)

    (42)

    NODE:=nops(eqode);NAE:=nops(eqae);

    NODE := 1

    NAE := 1

    (43)

    for XX from 1 to NODE do
    EQODE||XX:=lhs(eqode[XX])=rhs(eqode[XX])*ff:
    end do;

    EQODE1 := diff(y1(t), t) = z1(t)*(1/2+(1/2)*tanh(1000*t-1000))

    (44)

    for XX from 1 to NAE do
    EQAE||XX:=subs(pars,-mu*(diff(rhs(eqae[XX])-lhs(eqae[XX]),t))=rhs(eqae[XX])-lhs(eqae[XX]));
    end do;

    EQAE1 := .1*sin(z1(t))*(diff(z1(t), t))+.2*z1(t)*(diff(z1(t), t))+.1*(diff(z1(t), t))*(y1(t)+1)+.1*z1(t)*(diff(y1(t), t))+.1*(diff(y1(t), t)) = cos(z1(t))-z1(t)^2-z1(t)*(y1(t)+1)-y1(t)

    (45)

     

    Dvars1:={seq(diff(z||x(t),t)=D||x,x=1..NAE)};

    Dvars1 := {diff(z1(t), t) = D1}

    (46)

    Dvars2:={seq(rhs(Dvars1[x])=lhs(Dvars1[x]),x=1..NAE)};

    Dvars2 := {D1 = diff(z1(t), t)}

    (47)

    icsn:=seq(subs(y||x(0)=y||x(t),icodes[x]),x=1..NODE),seq(subs(z||x(0)=z||x(t),icaes[x]),x=1..NAE);

    icsn := y1(t) = 0., z1(t) = 0.

    (48)

    for j from 1 to NAE do

    EQAEX||j:=subs(Dvars1,eqode,icsn,Dvars2,lhs(EQAE||j))=rhs(EQAE||j);

    end do;

    EQAEX1 := .1*sin(0.)*(diff(z1(t), t))+.1*(diff(z1(t), t)) = cos(z1(t))-z1(t)^2-z1(t)*(y1(t)+1)-y1(t)

    (49)

    Sys:={seq(EQODE||x,x=1..NODE),seq(EQAEX||x,x=1..NAE),seq(icodes[x],x=1..NODE),seq(icaes[x],x=1..NAE)};

    Sys := {.1*(diff(z1(t), t)) = cos(z1(t))-z1(t)^2-z1(t)*(y1(t)+1)-y1(t), y1(0) = 0., z1(0) = 0., diff(y1(t), t) = z1(t)*(1/2+(1/2)*tanh(1000*t-1000))}

    (50)

    if Xexplicit=1 then

    sol:=dsolve(Sys,numeric):

    else

    sol:=dsolve(Sys,numeric,stiff=true,implicit=true):

    end if:

     

    for XX from 1 to NODE do
    a||XX:=odeplot(sol,[t,y||XX(t)],0..subs(pars,tf),color=red);
    end do:

    for XX from NODE+1 to NODE+NAE do
    a||XX:=odeplot(sol,[t,z||(XX-NODE)(t)],0..subs(pars,tf),color=blue);
    end do:

    display(seq(a||x,x=1..NODE+NAE),axes=boxed);

     

    End Example 3

     

    Example 4

    restart;

    with(plots):

    N:=11:h:=1/(N+1):

    for i from 2 to N+1 do eq1[i]:=diff(y||i(t),t)=(y||(i+1)(t)-2*y||i(t)+y||(i-1)(t))/h^2-y||i(t)*(1+z||i(t));od:

    for i from 2 to N+1 do eq2[i]:=0=(z||(i+1)(t)-2*z||i(t)+z||(i-1)(t))/h^2-(1-y||i(t)^2)*(exp(-z||i(t)));od:

    eq1[1]:=(3*y1(t)-4*y2(t)+y3(t))/(2*h)=0: eq1[N+2]:=y||(N+2)(t)-1=0:

    eq2[1]:=(3*z1(t)-4*z2(t)+1*z3(t))/(2*h)=0: eq2[N+2]:=z||(N+2)(t)=0:

    eq1[1]:=subs(y1(t)=z||(N+3)(t),eq1[1]):

    eq1[N+2]:=subs(y||(N+2)(t)=z||(N+4)(t),eq1[N+2]):

    eqode:=[seq(subs(y1(t)=z||(N+3)(t),y||(N+2)(t)=z||(N+4)(t),eq1[i]),i=2..N+1)]:

    eqae:=[eq1[1],eq1[N+2],seq(eq2[i],i=1..N+2)]:

    icodes:=[seq(y||j(0)=1,j=2..N+1)]:

    icaes:=[seq(z||j(0)=0,j=1..N+2),z||(N+3)(0)=1,z||(N+4)(0)=1]:

    pars:=[mu=0.00001,q=1000,tint=1,tf=2]:

    Xexplicit:=2:

    ff:=subs(pars,1/2+1/2*tanh(q*(t-tint))):

    NODE:=nops(eqode):NAE:=nops(eqae):

    for XX from 1 to NODE do

    EQODE||XX:=lhs(eqode[XX])=rhs(eqode[XX])*ff: end do:

    for XX from 1 to NAE do

    EQAE||XX:=subs(pars,-mu*(diff(rhs(eqae[XX])-lhs(eqae[XX]),t))=rhs(eqae[XX])-lhs(eqae[XX])); end do:

    Dvars1:={seq(diff(z||x(t),t)=D||x,x=1..NAE)}:

    Dvars2:={seq(rhs(Dvars1[x])=lhs(Dvars1[x]),x=1..NAE)}:

    icsn:=seq(subs(y||x(0)=y||x(t),icodes[x]),x=1..NODE),seq(subs(z||x(0)=z||x(t),icaes[x]),x=1..NAE):

    for j from 1 to NAE do

    EQAEX||j:=subs(Dvars1,eqode,icsn,Dvars2,lhs(EQAE||j))=rhs(EQAE||j):

    end do:

    Sys:={seq(EQODE||x,x=1..NODE),seq(EQAEX||x,x=1..NAE),seq(icodes[x],x=1..NODE),seq(icaes[x],x=1..NAE)}:

    if Xexplicit=1 then

    sol:=dsolve(Sys,numeric,maxfun=0):

    else

    sol:=dsolve(Sys,numeric,stiff=true,implicit=true,maxfun=0):

    end if:

     

    for XX from 1 to NODE do

    a||XX:=odeplot(sol,[t,y||(XX+1)(t)],1..subs(pars,tf),color=red): end do:

    for XX from NODE+1 to NODE+NAE do

    a||XX:=odeplot(sol,[t,z||(XX-NODE)(t)],1..subs(pars,tf),color=blue): end do:

    display(seq(a||x,x=1..NODE),a||(NODE+NAE-1),a||(NODE+NAE),axes=boxed);

     

    End of Example 4

     

    Sometimes the parameters of the switch function and perturbation need to be tuned to obtain propoer convergence. Below is Example 1 shown for several cases using the 'parameters' option in Maple's dsolve to compare how tuning parameters affects the solution

    restart:

    with(plots):

    eqode:=[diff(y1(t),t)=-y1(t)^2+z1(t)]: eqae:=[cos(y1(t))-z1(t)^0.5=0]:

    icodes:=[y1(0)=0.25]: icaes:=[z1(0)=0.8]:

    pars:=[tf=5]:

    Xexplicit:=2;

    Xexplicit := 2

    (51)

    ff:=subs(pars,1/2+1/2*tanh(q*(t-tint))):

    NODE:=nops(eqode):NAE:=nops(eqae):

    for XX from 1 to NODE do
    EQODE||XX:=lhs(eqode[XX])=rhs(eqode[XX])*ff:
    end do:

    for XX from 1 to NAE do
    EQAE||XX:=subs(pars,-mu*(diff(rhs(eqae[XX])-lhs(eqae[XX]),t))=rhs(eqae[XX])-lhs(eqae[XX]));
    end do:

     

    Dvars1:={seq(diff(z||x(t),t)=D||x,x=1..NAE)}:

    Dvars2:={seq(rhs(Dvars1[x])=lhs(Dvars1[x]),x=1..NAE)}:

    icsn:=seq(subs(y||x(0)=y||x(t),icodes[x]),x=1..NODE),seq(subs(z||x(0)=z||x(t),icaes[x]),x=1..NAE):

    for j from 1 to NAE do

    EQAEX||j:=subs(Dvars1,eqode,icsn,Dvars2,lhs(EQAE||j))=rhs(EQAE||j):

    end do:

    Sys:={seq(EQODE||x,x=1..NODE),seq(EQAEX||x,x=1..NAE),seq(icodes[x],x=1..NODE),seq(icaes[x],x=1..NAE)}:

    if Xexplicit=1 then

    sol:=dsolve(Sys,numeric,'parameters'=[mu,q,tint],maxfun=0):

    else

    sol:=dsolve(Sys,numeric,'parameters'=[mu,q,tint],stiff=true,implicit=true):

    end if:

     

    sol('parameters'=[0.1,1000,1]):

    plot1:=odeplot(sol,[t-1,y1(t)],0..4,color=red):
    plot2:=odeplot(sol,[t-1,z1(t)],0..4,color=blue):

    sol('parameters'=[0.001,10,1]):

    plot3:=odeplot(sol,[t-1,y1(t)],0..4,color=red,linestyle=dot):
    plot4:=odeplot(sol,[t-1,z1(t)],0..4,color=blue,linestyle=dot):

    display(plot1,plot2,plot3,plot4,axes=boxed);

     

    In general, one has to decrease mu, and increase q and tint until convergence (example at t=3)

    sol('parameters'=[0.001,10,1]):sol(3+1);

    [t = 4., y1(t) = .738587929442734, z1(t) = .546472878850096]

    (52)

    sol('parameters'=[0.0001,100,10]):sol(3+10);

    [t = 13., y1(t) = .738684397167344, z1(t) = .546618936273638]

    (53)

    sol('parameters'=[0.00001,1000,20]):sol(3+20);

    [t = 23., y1(t) = .738694113087217, z1(t) = .546633473784526]

    (54)

     

    The results have converged to 4 digits after the decimal. Of course, absolute and relative tolerances of the solvers can be modified if needed

     

    Download SolvingDAEsinMaple.mws

    I'd like to pay attention to Maple packages for symmetric functions, posets, root systems, and finite Coxeter groups by John Stembridge. That was estimated as well qualified in Mathoverflow.

    Dr Lopez,

    Thank you for your very interesting webinar on eigenpairs today. I especially valued your "aside" example of a matrix for which eigenvectors are "relatively simple," but which also has a "relatively simple" vector for which the magnitude is preserved, although not necessarily the direction. Your reply at the end of the program was very insightful, and I think it opens an interesting are of exploration.

    Juan Tolosa Richard Stockton College, NJ

    First 70 71 72 73 74 75 76 Last Page 72 of 301