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
  • I’m excited to announce that we’ve just launched MapleSim 2024.

    The new release has tools that are designed to drive innovation, and overall save you time when creating and developing simulations.

    At Maplesoft we are looking to continually enhance our engineering software with new features based on customer feedback, and I’m pleased to share some of the fruits of that labor, and thank the developers, product management team, and  customers that contributed.

    The new offering  helps engineers to

    • Rapidly Tune Parameters
    • Explore Design Concepts
    • Expand Modeling Capabilities

    For example, the new Rerun panel allows you to significantly cut the time between simulations as you quickly apply different parameter values, initial conditions and even simulation settings between runs. It does this by skipping the formulation steps when there are no structural changes made to the model, and you can even see the plots and results of the different iterations side by side.
    You can see it in action in this short demo video.

    There is now support for the latest Modelica feature set, so you can import Modelica Libraries that make use of MSL 4.0.0 features, and adds a range of new modeling components to the standard MapleSim libraries (electrical, 1-D mechanical, signal block and more).

     

    MapleSim 2024 also includes more components in the Hydraulic library to support modeling of flow restrictions and adds a Scripting button to add and organize Maple worksheets.

    We’ve also applied a whole series of updates to our MapleSim add-on products, including:

    • The MapleSim Web Handling Library has new tools for modeling heavier webs, winding of multiple rolls on the same drum, and adding a Switching Nip Roller to swing the web contact points between rollers.
    • The MapleSim Connector for FMI can now import and export FMI 3.0.
    • The MapleSim CAD Toolbox supports recent software releases from NX™, SOLIDWORKS®, Solid Edge®, Parasolid®, and other CAD tools.
    • The MapleSim Heat Transfer Library has gained a new T-junction component for the Water subpackage to improve flow/pressure-drop calculations for systems with branches.

    We have an upcoming webinar for you to see the new 2024 features in MapleSim Web Handling Library – you can sign up to register here.

    You can find out more about the other new features at the MapleSim What’s New web page, and as always, we are happy to hear your comments and product suggestions.

    And if you are new to MapleSim and would like to try building and running a model yourself, you can request a free trial, or contact Maplesoft sales team with any questions.

    Consider the equation  (2^x)*(27^(1/x)) = 24  for which we need to find the exact values ​​of its real roots. This is not difficult to solve by hand if you first take the logarithm of this equation to any base, after which the problem is reduced to solving a quadratic equation. But the  solve  command fails to solve this equation and returns the result in RootOf form. The problem is solved if we first ask Maple to take the logarithm of the equation. I wonder if the latest versions of Maple also do not directly address the problem?

    restart;
    Eq:=2^x*27^(1/x)=24:
    solve(Eq, x, explicit);
    
    map(ln, Eq); # Taking the logarithm of the equation
    solve(%, x);
    simplify({%}); # The final result
    

                      

     

    We've just launched Maple Flow 2024!

    You're in the driving seat with Maple Flow - each new feature has a straight-line connection to a user-driven demand to work faster and more efficiently.

    Head on over here for a rundown of everything that's new, but I thought I'd share my personal highlights here.

    If your result contains a large vector or matrix, you can now scroll to see more data. You can also change the size of the matrix to view more or fewer rows and columns.

    You can resize rows and columns if they're too large or small, and selectively enable row and column headers.

    If the vector or matrix in your result contains a unit, you can now rescale units with the Context Panel (for the entire matrix) or inline (for individual entries).

    A few releases ago, we introduced the Variables palette to help you keep track of all the user-defined parameters at point of the grid cursor.

    You can now insert variables into the worksheet from the Variables palette. Just double-click on the appropriate name.

    Maple Flow already features command completion - just type the first few letters of a command, and a list of potential completions appears. Just pick the completion you need with a quick tap of the Tab key.

    We've supercharged this feature to give potential arguments for many popular functions. Type a function name followed by an opening bracket, and a list appears.

    In case you've missed it, the argument completion list also features (when they make sense) user-defined variables.

    You can now link to different parts of the same worksheet. This can be used to create a table of contents that lets you jump to different parts of larger worksheets.

    This page lists everything that's new in the current release, and all the prior releases. You might notice that we have three releases a year, each featuring many user-requested items. Let me know what you want to see next - you might not have to wait that long!

    In our recent project, we're diving deep into understanding the SIR model—a fundamental framework in epidemiology that helps us analyze how diseases spread through populations. The SIR model categorizes individuals into three groups: Susceptible (S), Infected (I), and Recovered (R). By tracking how people move through these categories, we can predict disease dynamics and evaluate interventions.

    Key Points of the SIR Model:

    • Susceptible (S): Individuals who can catch the disease.
    • Infected (I): Those currently infected and capable of spreading the disease.
    • Recovered (R): Individuals who have recovered and developed immunity.

    Vaccination Impact: One of the critical interventions in disease control is vaccination, which moves individuals directly from the susceptible to the recovered group. This simple action reduces the number of people at risk, thereby lowering the overall spread of the disease.

    We're experimenting with a simple model to understand how different vaccination rates can significantly alter the dynamics of an outbreak. By simulating scenarios with varying vaccination coverage, students can observe how herd immunity plays a crucial role in controlling diseases. Our goal is to make these abstract concepts clear and relatable through practical modeling exercises.


     

    In this exercise, we are going back to the simple SIR model, without births or deaths, to look at the effect of vaccination. The aim of this activity is to represent vaccination in a very simple way - we are assuming it already happened before we run our model! By changing the initial conditions, we can prepare the population so that it has received a certain coverage of vaccination.

    We are starting with the transmission and recovery parameters  b = .4/daysand c = .1/days . To incorporate immunity from vaccination in the model, we assume that a proportion p of the total population starts in the recovered compartment, representing the vaccine coverage and assuming the vaccine is perfectly effective. Again, we assume the epidemic starts with a single infected case introduced into the population.​
    We are going to model this scenario for a duration of 2 years, assuming that the vaccine coverage is 50%, and plot the prevalence in each compartment over time.

     

    restart
    with(plots)

    b := .4; c := .1; n := 10^6; p := .5

    deS := diff(S(t), t) = -b*S(t)*I0(t); deI := diff(I0(t), t) = b*S(t)*I0(t)-c*I0(t); deR := diff(R(t), t) = c*I0(t)

    diff(R(t), t) = .1*I0(t)

    (1)

    F := dsolve([deS, deI, deR, S(0) = 1-p, I0(0) = 1/n, R(0) = p], [S(t), I0(t), R(t)], numeric, method = rkf45, maxfun = 100000)

    odeplot(F, [[t, S(t)], [t, I0(t)], [t, R(t)]], t = 0 .. 730, colour = [blue, red, green], legend = ["S(t)", "I0(t)", "R(t)"], labels = ["Time (days)", "  Proportion\nof Population "], title = "SIR Model with vaccine coverage 50 %", size = [500, 300])

     

    F(100)

    [t = 100., S(t) = HFloat(0.46146837378273076), I0(t) = HFloat(0.018483974421123688), R(t) = HFloat(0.5200486517961457)]

    (2)

    eval(S(:-t), F(100))

    HFloat(0.46146837378273076)

    (3)

    Reff := proc (s) options operator, arrow; b*(eval(S(:-t), F(s)))/(c*n) end proc; Reff(100)

    HFloat(1.845873495130923e-6)

    (4)

    plot(Reff, 0 .. 730, size = [500, 300])

     

    Increasing the vaccine coverage to 75%

    NULL

    restart
    with(plots)

    b := .4; c := .1; n := 10^6; p := .75

    deS := diff(S(t), t) = -b*S(t)*I0(t); deI := diff(I0(t), t) = b*S(t)*I0(t)-c*I0(t); deR := diff(R(t), t) = c*I0(t)

    diff(R(t), t) = .1*I0(t)

    (5)

    NULL

    F1 := dsolve([deS, deI, deR, S(0) = 1-p, I0(0) = 1/n, R(0) = p], [S(t), I0(t), R(t)], numeric, method = rkf45, maxfun = 100000)

    odeplot(F1, [[t, S(t)], [t, I0(t)], [t, R(t)]], t = 0 .. 730, colour = [blue, red, green], legend = ["S(t)", "I0(t)", "R(t)"], labels = ["Time (days)", "  Proportion\nof Population "], title = "SIR Model with vaccine coverage 75%", size = [500, 300])

     

    F(1100)

    eval(S(:-t), F1(100))

    HFloat(0.249990000844159)

    (6)

    Reff := proc (s) options operator, arrow; b*(eval(S(:-t), F1(s)))/(c*n) end proc; Reff(100)

    HFloat(9.99960003376636e-7)

    (7)

    plot(Reff, 0 .. 730, size = [500, 300])

     

    Does everyone in the population need to be vaccinated in order to prevent an epidemic?What do you observe if you model the infection dynamics with different values for p?

    No, not everyone in the population needs to be vaccinated in order to prevent an epidemic . In this scenario, if p equals 0.75 or higher, no epidemic occurs - 75 % is the critical vaccination/herd immunity threshold . Remember,, herd immunity describes the phenomenon in which there is sufficient immunity in a population to interrupt transmission . Because of this, not everyone needs to be vaccinated to prevent an outbreak .

    What proportion of the population needs to be vaccinated in order to prevent an epidemic if b = .4and c = .2/days? What if b = .6 and "c=0.1 days^(-1)?"

    In the context of the SIR model, the critical proportion of the population that needs to be vaccinated in order to prevent an epidemic is often referred to as the "herd immunity threshold" or "critical vaccination coverage."

    • 

    Scenario 1: b = .4and c = .2/days

    ``

    restart
    with(plots)

    b := .4; c := .2; n := 10^6; p := .5``

    deS := diff(S(t), t) = -b*S(t)*I0(t); deI := diff(I0(t), t) = b*S(t)*I0(t)-c*I0(t); deR := diff(R(t), t) = c*I0(t)

    diff(R(t), t) = .2*I0(t)

    (8)

    F1 := dsolve([deS, deI, deR, S(0) = 1-p, I0(0) = 1/n, R(0) = p], [S(t), I0(t), R(t)], numeric, method = rkf45, maxfun = 100000)

    odeplot(F1, [[t, S(t)], [t, I0(t)], [t, R(t)]], t = 0 .. 730, colour = [blue, red, green], legend = ["S(t)", "I0(t)", "R(t)"], labels = ["Time (days)", "  Proportion\nof Population "], title = "SIR Model with vaccine coverage 50 %", size = [500, 300])

     


    The required vaccination coverage is around 50% .

    • 

    Scenario 1: b = .6and c = .1/days

    restart
    with(plots)

    b := .6; c := .1; n := 10^6; p := .83NULL

    deS := diff(S(t), t) = -b*S(t)*I0(t); deI := diff(I0(t), t) = b*S(t)*I0(t)-c*I0(t); deR := diff(R(t), t) = c*I0(t)

    diff(R(t), t) = .1*I0(t)

    (9)

    NULL

    F1 := dsolve([deS, deI, deR, S(0) = 1-p, I0(0) = 1/n, R(0) = p], [S(t), I0(t), R(t)], numeric, method = rkf45, maxfun = 100000)

    odeplot(F1, [[t, S(t)], [t, I0(t)], [t, R(t)]], t = 0 .. 730, colour = [blue, red, green], legend = ["S(t)", "I0(t)", "R(t)"], labels = ["Time (days)", "  Proportion\nof Population "], title = "SIR Model with vaccine coverage 83% ", size = [500, 300])

     

    "The required vaccination coverage is around 83 `%` ."


    Download SIR_simple_vaccination_example.mw

    It can happen when an operation is interrupted by  that Maple does not return to  and still shows .

    This can give the false impression that the Maple server in charge of the evaluation did not get the message to stop whatever it was doing.

    By giving Maple an impossible task to solve analytically

    f1 := x1 - x1*sin(x1 + 5*x2) - x2*cos(5*x1 - x2);
    f2 := x2 - x2*sin(5*x1 - 3*x2) + x1*cos(3*x1 + 5*x2);
    solve({f1, f2});

    I have noticed in the Windows Task Manager that freeing allocated memory can take much longer than one might think.

    In one case it took 30 minutes to free 24 Gb of total allocated memory (21 Gb of it in RAM/physical memory). In this case the interrupt button became active (turned from grey to red ) two times and memory continued piling up  again.

    Lessons learned for me:

    • The task manager is not only a valuable indicator for task activity but also for the interruption/memory freeing process.
    • Before killing a whole Maple session and potentially losing the last state of a worksheet it can pay off to wait and repeatedly interrupt an operation.

     

    Suggestion: When the maple server gets an interrupt request, it could report to the GUI that it is in an interruption state and is no longer evaluating input. For example changing the message in the status bar from Evaluating... to Interrupting...


     

    When a derivative can be written as a function of the independent variable only for example

    y'=f(x)

    y''=f(x)

    y'''(x)=f(x)

    etc.

     

    We call that a directly integrable equation.

     

    Example 1:

     

    Find the general solution for the following directly integrable equation

    diff(y, x) = 6*x^2+4*y(1) and 6*x^2+4*y(1) = 0

    That means

    int(6*x^2+4, x)

    y = 2*x^3+c+4*x", where" c is an arbitary solution

    ``

     

     
    equation1 := diff(y(x), x) = 6*x^2+4

    diff(y(x), x) = 6*x^2+4

    (1)

    NULL

    NULL

    sol1 := dsolve(equation1, y(x))

    y(x) = 2*x^3+c__1+4*x

    (2)

    And  if we have the initial condition
    y(1) = 0
    particular_sol1 := dsolve({equation1, y(1) = 0}, y(x))

    y(x) = 2*x^3+4*x-6

    (3)

    "(->)"

     

     

     

     

     

    "Example 2:"NULL

    NULL

    "  Find the particular solution for the following equation with condition"

     

    x^2*(diff(y(x), x)) = -1

    y(1)=3

    So we will need to get the y' by itself

    int(-1/x^2, x)

    so,

    y = 1/x+c , where c is an arbitary constant

    And this is our general solution. Now we plug in the initial condition when x = 1, y = 3.

     

    That means c = 2.

     

    Thus, the particular solution is

     

    y = 1/x+2``

    eq := x^2*(diff(y(x), x)) = -1

    x^2*(diff(y(x), x)) = -1

    (4)

    NULL

    NULL

    sol := dsolve(eq, y(x))

    y(x) = 1/x+c__1

    (5)

    NULL

    particular_sol := dsolve({eq, y(1) = 3}, y(x))

    y(x) = 1/x+2

    (6)

    NULL

    NULL

    plot(1/x+2, x = -20 .. 20, color = "Red", axes = normal, legend = [typeset(1/x+2)])

     

    NULL

    NULL

    NULL

    NULL

    " Example 3:"

     

    " Find the particular solution for the following equation with condition"

     

    diff(y, t, t) = cost, (D(y))(0) = 0, y(0) = 1

    eq1 := diff(y(t), t, t) = cos(t)

    diff(diff(y(t), t), t) = cos(t)

    (7)

    particular_sol := dsolve({eq1, y(0) = 1, (D(y))(0) = 0}, y(t))

    y(t) = -cos(t)+2

    (8)

    "(->)"

     

     

     

     

    NULL


     

    Download integral.mw

     

    Hello,

    Attached I am sending several procedures for curves in 3D space. They were written without using Maple's built-in DiffGeo procedures and functions. As an example of their use, I made several animations - Maple worksheets are attached. I hope that maybe they will be useful to someone.

    Regards.

    ClsDGproc-Curves.zip

    restart;

    torus:= (x^2+y^2+z^2 + R^2-r^2)^2 - 4*R^2*(x^2+y^2):

    torus1:=eval(torus,[r=1,R=4]);

    (x^2+y^2+z^2+15)^2-64*x^2-64*y^2

    (1)

    plots:-implicitplot3d(torus1=0, x=-6..6,y=-6..6,z=-6..6, numpoints=5000, scaling=constrained, view=-2..2);

     

    sol1:=solve(torus1<0, [x,y,z]); # It should be easy for Maple

    [[x = 0, y < -3, -5 < y, z < (-y^2-8*y-15)^(1/2), -(-y^2-8*y-15)^(1/2) < z], [x = 0, y < 5, 3 < y, z < (-y^2+8*y-15)^(1/2), -(-y^2+8*y-15)^(1/2) < z]]

    (2)

    eval(torus1, [x=3,y=2,z=0]); # in the interior

    -48

    (3)

    eval(sol1, [x=3,y=2,z=0]);   # ???

    [[3 = 0, 2 < -3, -5 < 2, 0 < (-35)^(1/2), -(-35)^(1/2) < 0], [3 = 0, 2 < 5, 3 < 2, 0 < (-3)^(1/2), -(-3)^(1/2) < 0]]

    (4)

     

     

    Download Bug-solve_ineqs.mw

    This Maplesoft guest blog post is from Prof. Dr. Johannes Blümlein from Deutsches Elektronen-Synchrotron (DESY), one of the world’s leading particle accelerator centres used by thousands of researchers from around the world to advance our knowledge of the microcosm. Prof. Dr. Blümlein is a senior researcher in the Theory Group at DESY, where he and his team make significant use of Maple in their investigations of high energy physics, as do other groups working in Quantum Field Theory. In addition, he has been involved in EU programs that give PhD students opportunities to develop their Maple programming skills to support their own research and even expand Maple’s support for theoretical physics.


     

    The use of Maple in solving frontier problems in theoretical high energy physics

    For several decades, progress in theoretical high energy physics relies on the use of efficient computer-algebra calculations. This applies both to so-called perturbative calculations, but also to non-perturbative computations in lattice field theory. In the former case, large classes of Feynman diagrams are calculated analytically and are represented in terms of classes of special functions. In early approaches started during the 1960s, packages like Reduce [1] and Schoonship [2] were used. In the late 1980s FORM [3] followed and later on more general packages like Maple and Mathematica became more and more important in the solution of these problems. Various of these problems are related to data amounts in computer-algebra of O(Tbyte) and computation times of several CPU years currently, cf. [4].

    Initially one has to deal with huge amounts of integrals. An overwhelming part of them is related by Gauss’ divergence theorem up to a much smaller set of the so-called master integrals (MIs). One performs first the reduction to the MIs which are special multiple integrals. No general analytic integration procedures for these integrals exist. There are, however, several specific function spaces, which span these integrals. These are harmonic polylogarithms, generalized harmonic polylogarithms, root-valued iterated integrals and others. For physics problems having solutions in these function spaces codes were designed to compute the corresponding integrals. For generalized harmonic polylogarithms there is a Maple code HyperInt [5] and other codes [6], which have been applied in the solution of several large problems requiring storage of up to 30 Gbyte and running times of several days. In the systematic calculation of special numbers occurring in quantum field theory such as the so-called β-functions and anomalous dimensions to higher loop order, e.g. 7–loop order in Φ4 theory, the Maple package HyperLogProcedures [7] has been designed. Here the largest problems solved require storage of O(1 Tbyte) and run times of up to 8 months. Both these packages are available in Maple only.

    A very central method to evaluate master integrals is the method of ordinary differential equations. In the case of first-order differential operators leading up to root-valued iterative integrals their solution is implemented in Maple in [8] taking advantage of the very efficient differential equation solvers provided by Maple. Furthermore, the Maple methods to deal with generating functions as e.g. gfun, has been most useful here. For non-first order factorizing differential equation systems one first would like to factorize the corresponding differential operators [9]. Here the most efficient algorithms are implemented in Maple only. A rather wide class of solutions is related to 2nd order differential equations with more than three singularities. Also here Maple is the only software package which provides to compute the so-called 2F1 solutions, cf. [10], which play a central role in many massive 3-loop calculations

    The Maple-package is intensely used also in other branches of particle physics, such as in the computation of next-to-next-to leading order jet cross sections at the Large Hadron Collider (LHC) with the package NNLOJET and double-parton distribution functions. NNLOJET uses Maple extensively to build the numerical code. There are several routines to first build the driver with automatic links to the matrix elements and subtraction terms, generating all of the partonic subprocesses with the correct factors. To build the antenna subtraction terms, a meta-language has been developed that is read by Maple and converted into calls to numerical routines for the momentum mappings, calls to antenna and to routines with experimental cuts and plotting routines, cf. [11].

    In lattice gauge calculations there is a wide use of Maple too. An important example concerns the perturbative predictions in the renormalization of different quantities. Within different European training networks, PhD students out of theoretical high energy physics and mathematics took the opportunity to take internships at Maplesoft for several months to work on parts of the Maple package and to improve their programming skills. In some cases also new software solutions could be obtained. Here Maplesoft acted as industrial partner in these academic networks.

    References

    [1] A.C. Hearn, Applications of Symbol Manipulation in Theoretical Physics, Commun. ACM 14 No. 8, 1971.

    [2] M. Veltman, Schoonship (1963), a program for symbolic handling, documentation, 1991, edited by D.N. Williams.

    [3] J.A.M. Vermaseren, New features of FORM, math-ph/0010025.

    [4] J. Blümlein and C. Schneider, Analytic computing methods for precision calculations in quantum field theory, Int. J. Mod. Phys. A 33 (2018) no.17, 1830015 [arXiv:1809.02889 [hep-ph]].

    [5] E. Panzer, Algorithms for the symbolic integration of hyperlogarithms with applications to Feynman integrals, Comput. Phys. Commun. 188 (2015) 148–166 [arXiv:1403.3385 [hep-th]].

    [6] J. Ablinger, J. Blümlein, C .Raab, C. Schneider and F. Wissbrock, Calculating Massive 3-loop Graphs for Operator Matrix Elements by the Method of Hyperlogarithms, Nucl. Phys. 885 (2014) 409-447 [arXiv:1403.1137 [hep-ph]].

    [7] O. Schnetz, φ4 theory at seven loops, Phys. Rev. D 107 (2023) no.3, 036002 [arXiv: 2212.03663 [hep-th]].

    [8] J. Ablinger, J. Blümlein, C. G. Raab and C. Schneider, Iterated Binomial Sums and their Associated Iterated Integrals, J. Math. Phys. 55 (2014) 112301 [arXiv:1407.1822 [hep-th]].

    [9] M. van Hoeij, Factorization of Differential Operators with Rational Functions Coefficients, Journal of Symbolic Computation, 24 (1997) 537–561.

    [10] J. Ablinger, J. Blümlein, A. De Freitas, M. van Hoeij, E. Imamoglu, C. G. Raab, C. S. Radu and C. Schneider, Iterated Elliptic and Hypergeometric Integrals for Feynman Diagrams, J. Math. Phys. 59 (2018) no.6, 062305 [arXiv:1706.01299 [hep-th]].

    [11] A. Gehrmann-De Ridder, T. Gehrmann, E.W.N. Glover, A. Huss and T.A. Morgan, Precise QCD predictions for the production of a Z boson in association with a hadronic jet, Phys. Rev. Lett. 117 (2016) no.2, 022001 [arXiv:1507.02850 [hep-ph]].

         Happy Easter to all those who celebrate! One common tradition this time of year is decorating Easter eggs. So, we’ve decided to take this opportunity to create some egg-related math content in Maple Learn. This year, a blog post by Tony Finch inspired us to create a walkthrough exploring the four-point egg. The four-point egg is a method to construct an egg-shaped graph using just a compass and a ruler, or in this case, Maple Learn. Here's the final product: 

         The Maple Learn document, found here, walks through the steps. In general, each part of the egg is an arc corresponding to part of a circle centred around one of the points generated in this construction. 

         For instance, starting with the unit circle and the three red points in the image below, the blue circle is centred at the bottom point such that it intersects with the top of the unit circle, at (0,1). The perpendicular lines were constructed using the three red points, such that they intersect at the bottom point and pass through opposite side points, either (-1,0) or (1,0). Then, the base of the egg is constructed by tracing an arc along the bottom of the blue circle, between the perpendicular lines, shown in red below.

     

         Check out the rest of the steps in the Maple Learn Document. Also, be sure to check out other egg-related Maple Learn documents including John May’s Egg Formulas, illustrating other ways to represent egg-shaped curves with mathematics, and Paige Stone’s Easter Egg Art, to design your own Easter egg in Maple Learn. So, if you’ve had your fill of chocolate eggs, consider exploring some egg-related geometry - Happy Easter!  

    Here are the Maple sources for the paper "Maximum Gap among Integers Having a Common Divisor with an Odd semi-prime" published on Journal of Advances in Mathematics and Computer Science 39 (10):51-61, at :https://doi.org/10.9734/jamcs/2024/v39i101934.

     

    gap := proc(a, b) return abs(a - b) - 1; end proc

    HostsNdivisors := proc(N)

    local i, j, g, d, L, s, t, m, p, q, P, Q, np, nq;

    m := floor(1/2*N - 1/2);

    L := evalf(sqrt(N));

    P := Array();

    Q := Array();

    s := 1; t := 1;

    for i from 3 to m do

       d := gcd(i, N);

        if 1 < d and d < L then P(s) := i; s++;

        elif L < d then Q(t) := i; t++; end if;

    end do;

      np := s - 1;

      nq := t - 1;

     for i to np do printf("%3d,", P(i)); end do;

      printf("\n");

      for i to nq do printf("%3d,", Q(i)); end do;

      printf("\n gaps: \n");

      for i to np do

         for j to nq do

          p := P(i); q := Q(j);

          g := gap(p, q);

          printf("%4d,", g);

      end do;

        printf("\n");

    end do;

    end proc

     

    HostOfpq := proc(p, q)

    local alpha, s, t, g, r, S, T, i, j;

       S := 1/2*q - 1/2;

       T := 1/2*p - 1/2;

       alpha := floor(q/p);

        r := q - alpha*p;

       for s to S do

         for t to T do

           g := abs((t*alpha - s)*p + t*r) - 1;

            printf("%4d,", g);

          end do;

         printf("\n");

     end do;

    end proc

    MapleSource.pdf

    MapleSource.mw

     

    To celebrate this day of mathematics, I want to share my favourite equation involving Pi, the Bailey–Borwein–Plouffe (BBP) formula:

    This is my favourite for a number of reasons. Firstly, Simon Plouffe and the late Peter Borwein (two of the authors that this formula is named after) are Canadian! While I personally have nothing to do with this formula, the fact that fellow Canadians contributed to such an elegant equation is something that I like to brag about.

    Secondly, I find it fascinating how Plouffe first discovered this formula using a computer program. It can often be debated whether mathematics is discovered or invented, but there’s no doubt here since Plouffe found this formula by doing an extensive search with the PSLQ integer relation algorithm (interfaced with Maple). This is an example of how, with ingenuity and creativity, one can effectively use algorithms and programs as powerful tools to obtain mathematical results.

    And finally (most importantly), with some clever rearranging, it can be used to compute arbitrary digits of Pi!

    Digit 2024 is 8
    Digit 31415 is 5
    Digit 123456 is 4
    Digit 314159 is also 4
    Digit 355556 is… F?

    That last digit might look strange… and that’s because they’re all in hexadecimal (base-16, where A-F represent 10-15). As it turns out, this type of formula only exists for Pi in bases that are powers of 2. Nevertheless, with the help of a Maple script and an implementation of the BBP formula by Carl Love, you can check out this Learn document to calculate some arbitrary digits of Pi in base-16 and learn a little bit about how it works.

    After further developments, this formula led to project PiHex, a combined effort to calculate as many digits of Pi in binary as possible; it turns out that the quadrillionth bit of Pi is zero! This also led to a class of BBP-type formulas that can calculate the digits of other constants like (log2)*(π^2) and (log2)^5.

    Part of what makes this formula so interesting is human curiosity: it’s fun to know these random digits. Another part is what makes mathematics so beautiful: you never know what discoveries this might lead to in the future. Now if you’ll excuse me, I have a slice of lemon meringue pie with my name on it 😋

     

    References
    BBP Formula (Wikipedia)
    A Compendium of BBP-Type Formulas
    The BBP Algorithm for Pi

    This post is to help anyone who is just as frustrated about typesetting in plots
    as I was before I solved my problem. (Note: the technique works with the version

    2020 and newer and may work with earlier versions.)

     

    Why? Because there is no obvious help in Maple showing what works. (And if folks

    can improve what I have posted,  please do so. At least when someone executes

    a search for this type of problem, they might see the best approach.)

    Goal: typeset a name of a function in any text of a plot.

    Approach: According to help, '?plot, typesetting',  one should use the
    option(procedure?) typeset.   For example:

    "restart;   plots:-setoptions(size = [400, 200]):    `f__1`(x) :=cos(x)*(e)^(-(x^(2))/(4)) :    p1 := plot(`f__1`(x), x = -5..5,                         legend = typeset("function ", `f__1`(x) )  );"

     

     

    However, what I want in the legend is the expression "`f__1`(x),"not the evaluated

    expression. Entering the name with single quotes around the expression leads to this:

     

    p1 := plot(f__1(x), x = -5 .. 5, legend = typeset("function ", 'f__1(x)'))

     

     

    which is great, except that when I wish to redisplay the plot

     

    plots:-display(p1)

     

    the expression f__1(x)is evaluated.

     

    According to the code completion capability of Maple, the procedure

    "Typesetting:-Typeset" exists, and it does not evaluate the function:

     

    p1 := plot(f__1(x), x = -5 .. 5, legend = Typesetting:-Typeset(f__1(x))); plots:-display(p1)

     

     

    except there is no help regarding this procedure.

     

    It appears that the procedure operates only on one item.

     

    Solution: Hence, the ultimate solution for my problem is to still use the typeset 

    option, but Typeset any expression.

     

    p1 := plot(f__1(x), x = -5 .. 5, legend = typeset("function ", Typesetting:-Typeset(f__1(x)))); plots:-display(p1)

     

     

    Again, if you can contribute to this post, I would appreciate it.


     

    Download MaplePrimies_-_Typesetting_in_plots.mw



    (EDITED 2024/03/11  GMT 17H)

    In a recent Question@cq mentionned in its last reply "In fact, I wanted to do parameter sensitivity analysis and get the functional relationship between the [...] function and [parameters]. Later, i will study how the uncertainty of [the parameters] affects the [...] function".
    I did not keep exchanging further on with @cq, simply replying that I could provide it more help if needed.

    In a few words the initial problem was this one:

    • Let X_1 and X_2 two random variables and G the random variable defined by  G = 1 - (X_1 - 1)^2/9 - (X_2 - 1)^3/16.
       
    •  X_1 and X_2 are assumed to be gaussian random variables with respective mean and standard deviation equal to (theta_1, theta_3) and (theta_2, theta_4).
       
    • The four theta parameters are themselves assumed to be realizations of four mutually independent uniform random variables Theta_1, ..., Theta_4 whose parameters are constants.
       
    • Let QOI  (Quantity Of Interest) denote some scalar statistic of G (for instance its Mean, Variance, Skweness, ...).
      For instance, if QOI = Mean(G), then  QOI expresses itself as a function of the four parameters theta_1, ..., theta_4.
      The goal of @cq is to understand which of those parameters have the greatest influence on QOI.


    For a quick survey of Sensitivity Analysys (SA) and a presentation of some of the most common strategies see Wiki-Overview


    The simplest SA is the Local SA (LSA) we are all taught at school: having chosen some reference point P in the [theta_1, ..., theta_4] space the 1st order partial derivatives d[n] = diff(QOI, theta_n) expressed at point P give a "measure" (maybe after some normalization) of the sensitivity, at point P, of QOI regardibg each parameter theta_n.


    A more interesting situation occurs when the parameters can take values in a neighorhood of  P which is not infinitesimal, or more generally in some domain without reference to any specific point P.
    That is where Global SA (GSA) comes into the picture.
    While the notion of local variation at some point P is well established, GSA raises the fundamental question of how to define how to measure the variation of a function over an arbitrary domain?
    Let us take a very simple example while trying to answer this question "What is the variation of sin(x) over [0, 2*Pi]?"

    1. If we focus on the global trend of sin(x)  mean there is no variation at all.
    2. If we consider peak-to-peak amplitude the variation is equal to 2.
    3. At last, if we consider L2 norm the variation is equal to Pi.
      (but the constant function x -> A/sqrt(2) has the same L2 norm but it is flat, and in some sense les fluctuating). 


    Statisticians are accustomed to use the concept of variance as a measure to quantify the dispersion of a random variable. At the end of the sixties  one of them, Ilya Meyerovich Sobol’,  introduced the notion of Variance-Based GSA as the key tool to define the global variation of a function. This notion naturally led to that of Sobol' indices as a measure of the sensitivity of a funcion regarding one of its parameters or, which most important, regarding any combination (on says interaction) of its parameters.

    The aim of this post is to show how Sobol' indices can be computed when the function under study has an analytic expression.
     

    The Sobol' analysis is based on an additive decomposition of this function in terms of 2^P mutually orthogonal fiunctions where P is the number of its random parameters.
    This decomposition and the ensuing integrations whose values will represent the Sobol' indices can be done analytically in some situation. When it is no longer the case specific numerical estimation methods have to be used;


    The attached file contains a quite generic procedure to compute exact Sobol' indices and total Sobol' indices for a function whose parameters have any arbitrary statistical distribution.
    Let's immediately put this into perspective by saying that these calculations are only possible if Maple is capable to find closed form expressions of some integrals, which is of course not always the case.

    A few examples are also provided, including the one corresponding to @cq's original question.
    At last two numerical estimation methods are presented.

    SOBOL.mw

     

         On International Women’s Day we celebrate the achievements of women around the world. One inspiring story of women in STEM is that of Sophie Germain (1776-1831), a French mathematician and physicist who made groundbreaking strides in elasticity theory and number theory. She learned mathematics from reading books in her father’s library, and while she was not permitted to study at the École Polytechnique, due to prejudice against her gender, she was able to obtain lecture notes and decided to submit work under the name Monsieur LeBlanc. Some prominent mathematicians at the time, including Joseph-Louis Lagrange and Carl Friedrich Gauss, with whom Germain wrote, recognized her intellect and were supportive of her work in mathematics. 

         Sophie Germain is remembered as a brilliant and determined trailblazer in mathematics. She was the first woman to win a prize from the Paris Academy of Sciences for her contributions in elasticity theory and was among the first to make significant contributions toward proving Fermat’s Last Theorem. Among her many accomplishments, one special case of Fermat’s Last Theorem she was able to prove is when the exponent takes the form of what is now known as a Sophie Germain prime: a prime, p, such that 2p+1 is also a prime. The associated prime, 2p+1, is called a safe prime. 

         To mark International Women’s Day, I’ve created a document exploring the Ulam spiral highlighting Sophie Germain primes and safe primes, as an adaptation of Lazar Paroski’s Ulam spiral document. The image below displays part of the Ulam spiral with Sophie Germain primes highlighted in red, safe primes highlighted in blue, primes that are both a Sophie Germain prime and safe prime highlighted in purple, and primes that are neither in grey. 

      

         The document also contains small explorations of these types of prime numbers. For instance, one interesting property of safe primes is that they must either be 5, 7 or take the form 12k-1 for some k≥1. This can be shown from the fact a safe prime q must equal 2p+1 for some prime, p (a Sophie Germain prime), by definition. Then, since q and p are prime, for q > 7 we can determine through contradiction that q ≡ 3 (mod 4) and q ≡ 5 (mod 6), to conclude q ≡ 11 (mod 12) ≡ -1 (mod 12). And so, q = 12k-1 for some k≥1. The Maple Learn document can be found here along with its Maple script. The document also includes a group where you can test whether some positive integer of your choice, n, is a Sophie Germain prime or a safe prime. Alternatively, given n, a button press will display the next Sophie Germain prime greater than n, using Maple’s NextSafePrime command in the number theory package.  

         In mathematics, there is no shortage of interesting rabbit holes to dive into; many of which are the result of past and present women in mathematics, like Sophie Germain, who have persevered despite their hardships. 

    1 2 3 4 5 6 7 Last Page 3 of 300