The Lunar New Year is approaching and 2024 is the Year of the Dragon! This inspired me to create a visualization approximating the dragon curve in Maple Learn, using Maple. 

The dragon curve, first described by physicist John Heighway, is a fractal that can be constructed by starting with a single edge and then continually performing the following iteration process:  

Starting at one endpoint of the curve, traverse the curve and build right triangles on alternating sides of each edge on the curve. Then, remove all the original edges to obtain the next iteration. 

This process continues indefinitely, so while we can’t draw the fractal perfectly, we can approximate it using a Lindenmayer system. In fact, Maple can do all the heavy lifting with the tools found in the Fractals package, which includes the LSystem subpackage to build your own Lindenmayer systems. The subpackage also contains different examples of fractals, including the dragon curve. Check out the Maple help pages here: 

Overview of the Fractals Package  

Overview of the Fractals:-LSystem Subpackage 

Using this subpackage, I created a Maple script (link) to generate a Maple Learn document (link) to visualize the earlier iterations of the approximated dragon curve. Here’s what iteration 11 looks like: 

You can also add copies of the dragon curve, displayed at different initial angles, to visualize how they can fit together. Here are four copies of the 13th iteration: 

Mathematics is full of beauty and fractals are no exception. Check out the LSystemExamples subpackage to see many more examples. 

Happy Lunar New Year! 

A ball on a turntable can move in circles instead of falling off the edge (provided the initial conditions are appropriate). The effect was demonstrated in a video and can be simulated with MapleSim. The amination below shows a simulation of a frictionless case (falling off the table) and the case with a coefficient of friction of one.

Also demonstrated in the video: Tilting the table leads to a sideward (not a downhill) movement of the ball.

The presenter of the video noted that in the untilted state, the ball eventually drifts off the table, attributing this to slippage. This drift is also observable in the animation above, where the ball starts moving in a spiral, whereas in a Maple simulation (below to the left), the ball follows a perfect circle. Only after optimizing contact and initial conditions, MapleSim produced a result (using the same parameters) that exhibits a similar circle, with a slight difference in angular orientation after completing two revolutions about the center of the circle.

Some observations on the MapleSim model:

• The results only slightly depend on the solvers. Numerical inaccuracies do not seem to be the reason for the difference in orientation. (Edit: see update below for the reason).
• The ball bounces up and down in the MapleSim simulation (0.0025 of the balls radius). The bouncing is caused by the fact that the initial position of the ball is above the elastic equilibrium position with respect to the table (the elastic contact makes the ball sink in a bit). Adding damping in the settings of the contact element attenuates this effect and reduces the drift.
• Drift is not observable anymore if in the contact element setting for "kmu" (smoothness coefficient of sliding friction) is set to larger values (above 10 in this example). This is an indication that sliding friction occurs during the simulation.
• Choosing the equilibrium position as initial condition for the ball does not prevent initial bouncing because MapleSim sets an initial velocity for the ball that is directed away from the table. I did not manage to enforce strictly zero velocity. Maybe someone can tell why that is and how to set MapleSim to start the simulation without vertical velocity. (Edit: see update below for the reason).
• Assuming an initial velocity towards the contact to cancel the initial vertical velocity set by MapleSim substantially reduces bouncing and increases the diameter of the circle. This finally leads to a diameter that matches the Maple simulation. Therefore the initial bouncing combined with slippage seems to dissipate energy which leads to smaller circles. Depending on the contact settings and initial conditions for vertical movement the diameter of the circle varied moderately by about 15%.

In summary, MapleSim can be parametrized to simulate an ideal case without slippage. From there it should be possible to tune contact parameters to closely reproduce experiments where parameters are often not well known in advance.  

Some thoughts for future enhancement of MapleSim:

• In the model presented here, a patterned ball would have been helpful to visualize the tumbling movement of the ball. Marking two opposite sides of the ball with colored smaller spheres is a fiddly exercise that does not look nice.
• A sensor component that calculates the energy of a moving rigid body would have helped identifying energy loss of the system. Ideally the component could have two ports for the rotational and translational energy components. I see professional interest in such calculations and not only educational value for toy experiments.
• A port for slippage on the contact elements would have helped understanding where slippage occurs. Where slippage is, there is wear and this is also of interest for industrial applications.

Turntable_Paradox.msim (contains parameter sets for the above observations)

February 2nd is Groundhog Day in North America. A day when we look to small marmots to prognosticate the weather. If the groundhog sees its shadow when it emerges from its burrow, then it predicts 6 more weeks of winter, and if not, then spring begins today! Unfortunately there are many official weather predicting groundhogs. Fortunately, the excellent website Groundhog-Day.com tracks each of their predictions. Unfortunately, it doesn't tell you which groundhog to trust. Fortunately, it has an API and we can take the data and map it in Maple:

This map assumes that each groundhog's prediction is valid at it's exact geographic coordinates, but that it's predictive powers fall off in inverse proportion to the distance away.  So, exactly halfway between a groundhog predicting early spring, and one predicting 6 more weeks of winter, we expect 3 more weeks of winter.  I handle that with Maple's Interpolation:-InverseDistanceWeightedInterpolation command with a radius of 1500 miles.  I plot a contourplot of that interpolating function, and then display it with the world map in DataSets:-Builtin:-WorldMap to generate the image above.

All the code to do that can be found in the following worksheet which also using the URL packaget to fetch the most recent groundhog data possible from the website.

Groundhog-Map.mw

I've commented out a few lines that you might use to explore other possible maps.  You can filter to file to just include real living groundhogs and not all the other precitions (some from puppets, some from other animals) if you find that more trust worthy. You might also prefer to change the interpolation command, one of my collegues suggests that Interpolation:-NaturalNeighborInterpolation might be a better choice.

[Right-click on image and open in new tab to see larger]

To Scan Math with the Maple Calculator and show solution steps in Maple:

1. first scan some math with the calculator

2. Maple calculator immediately shows the solution if that is what you are looking for:

3. Calculator gives options to show the solution steps in the calculator itself ( footprint button in top-right) 

4. Or to upload the math to the MapleCloud (cloud icon with up arrow)

5. Once the math is uploaded, MapleCloud can be loaded on a desktop computer and the file opened from your account's Maple Calculator group of files:

6. Again, the solution and some more details are visible on Maple Cloud:

7. To open this math in Maple, click the blue button to Download the file.

The downloaded file can then be loaded in Maple:

8. The Maple commands to solve this math are shown, and the result. 

To show steps in Maple at this point, convert the math to inert form, then run the Student:-Calculus1:-ShowSolution() command on it:

Download MapleCalculatorMathCloudUpload.mw

There appears to be a bug with Maple 2023.2 which will remove units from physical constants.

AddConstant(Solar_equatorial_radius, symbol = r[e,Sol], value = 696342., uncertainty = 65., units = km) :

AddConstant(Solar_flattening, symbol = f[Sol], value = 0.000009) :

AddConstant(Solar_polar_radius, symbol = r[p,Sol], derive = -r[e,Sol]*(f[Sol] - 1)) :

AddConstant(Solar_nonradius, symbol = x[Sol], derive = f[Sol]*r[e,Sol]) :

In the example here, "GetUnit(Constant(r[p,Sol]))" will return "1" and "GetUnit(Constant(x[Sol]))" will return "m".

The only workaround is to not save dimensionless quantities as physical constants or otherwise at least some constants derived from a dimensionless constant will also be dimensionless.

Adding "units = 1" to "AddConstant" does not help.

The original app center was excellent, they were .. "finessed" so to speak.  Then they changed it (like so many websites today most are clunky and they're all terrible, they don't work well, and I'm not just singling out maplesoft)  Most websites that are finely tuned are most likely based on old school programming, just recall the old days of mapleprimes - the forum was nicely done - but that's a debate for later. Now to the nitty gritty.

I did a search at maplesoft app center of a particular author, and I know he has at least 20 applications.  Only 3 came up, and the last app entry said more apps by this author but it kept disappearing and I kept getting thrown back to application #1 all the meanwhile looking like I'm scrolling further down in the list. 

Maybe this is browser issue?  I'm using firefox.  Maybe Microsoft Edge works better? 

I'm an engineer and when showing results of calculations, some values will display as fractions, and I would prefer that instead floating numbers are displayed.  Also, there is kind of a quirk where if the multiplier of a unit is 1, the result displays as a unit only. I would prefer to see 1*A rather than A.

I wrote this simple proc to convert a value with or without attached unit to a floating point number if it is a fraction or if it has a unit and the coefficient is 1.

Let me know if there is a more elegant way to do this or you have any suggestions or questions.

   unrationalize := proc(x)
        local 
            returnval,
            localcoeff,
            localunit
        ;
        description
            "Converts a fractional number to float",
            "Units are supported"
        ;
        if type(x, fraction) or type(x, with_unit(fraction)) then 
            returnval := evalf(x)
        elif type(x, with_unit(1, 'localcoeff', 'localunit')) then
            returnval :=  evalf(localcoeff)*localunit
        else 
            returnval := x;
        end if;
        return returnval;
    end  proc;
# Testing the proc
list1 := [1/2, 1/2*Unit(('A')/('V')), 1, Unit('A')];
listDescription := ["Fraction", "Fraction with Unit", "Unity", "Unity with Unit"];
for i, myValue in list1 do
    [listDescription[i], "evalf:", evalf(myValue), "unrationalize:", unrationalize(myValue)];
end do;

Curve sketching is an important skill for all calculus students to learn. In an era where technology is increasingly relied upon to perform mathematical computations and representations, maintaining fundamental skills such as curve sketching is more important than ever.

The new “Curve Sketching” collection is now available on Maple Learn. This collection provides background information on the process of curve sketching and opportunities to put this knowledge into practice. By starting with the “Curve Sketching Guide” and “Relationships Between Derivatives” documents, students are exposed to observational and computational strategies for drawing a function and its 1st and 2nd derivatives.

After looking through these documents, students can begin to practice sketching by observing and interpreting graphical properties with the “Sketch Derivative From Function Graph”, “Sketch Second Derivative From Function Graph”, and “Sketch Function From Derivative Graphs” activities:

Once a student has mastered extracting sketching information by graphical observation, they are ready for the next step: extracting information from a function’s definition. At this point, the student is ready to try sketching from a blank canvas with the “Sketch Curve From Function Definition” activity:

This collection also has activities for students below the calculus level. For example, the “Curve Sketching Quadratics Activity”, can be completed using only factoring strategies:

Whether you are a quadratics rookie or a calculus pro, this collection has an interactive activity to challenge your knowledge. Have fun sketching!

My friend and colleague Nic Fillion and I are writing another book, this one on perturbation methods using backward error analysis (and Maple).  We have decided to make the supporting materials available by means of Jupyter notebooks with a Maple kernel (there are some Maple worksheets and workbooks already, but going forward we will use Jupyter).

The presentation style is meant to aid reproducibility, and to allow others to solve related problems by changing the scripts as needed.

The first one is up at 

https://github.com/rcorless/Perturbation-Methods-in-Maple

Comments very welcome.  This particular method is a bit advanced in theory (but it's very simple in practice, for weakly nonlinear oscillators).  I haven't coded for efficiency and there may be some improvements possible ("may" he says, sheesh).  Comments on that are also welcome.

-r

In the most recent issue of Maple Transactions, I published (with David Jeffrey, and with a student named Johan Joby) a paper that used Jupyter Notebook with a Maple kernel as the main vehicle.  Have a look, and let me know what you think.

Two-cycles in the infinite exponential tower

Happy new year everyone. Instead on making my usual new year resolutions, I though I would put up a request for some helpful  improvements here.

1) Please make search work properly. Allow searches to target specific sections. Posts, Questions, My own account, Other users. etc,

2) Add links to Maple help. Is there an online way of viewing Maple document? Or build Maple viewer into this site.

3) Add links to Maple Learn.

Download bugs-irrationals.mw

Two solstices occur on Earth every year, around June 21^st and December 21^st, often called the “June Solstice” and the “December Solstice” respectively. These solstices occur when the sun reaches its northernmost or southernmost point relative to the equator. During a solstice, the Northern Hemisphere will either experience the most sunlight of the year or the least sunlight of the year, while the Southern Hemisphere will experience the opposite phenomenon. The hemisphere with the most sunlight experiences a summer solstice, while the other hemisphere experiences a winter solstice.

Canada is located in the Northern Hemisphere and this Thursday, December 21^st, we will be experiencing a winter solstice. As the day with the least sunlight, this will be the shortest day of the year and consequently the longest night of the year.

Here in Canada, the sun will reach its minimum elevation during the winter solstice, and it will reach its maximum elevation during the Southern Hemisphere’s summer solstice on the same day. 

How high in the sky does the sun really get during these solstices? Check out our new Maple Learn document, Winter and Summer Solstice Sun Angles to find out. The answer depends on your latitude; for instance, with a latitude of approximately 43.51°, the document helps us find that the maximum midday elevation of the sun, which occurs during a summer solstice, will be 69.99°.

But how is the latitude of a location determined in the first place? See Maple Learn’s Calculating Latitude document to find out how the star Polaris, the center of the Earth, and the equator are all connected to latitude.

Latitude is one of two geographical coordinates that are paired together to specify a position on Earth, the other being longitude. See our Calculating Longitude document to explore how you can use your local time to approximate your longitude.

Armed with these coordinates, you can describe your position on the planet and find any number of interesting facts, such as your solstice sun angles from earlier, the time for sunrise and sunset, and the position of the Moon.

Happy Winter Solstice!

Why this post
This work was intended to be a simple reply to a question asked a few days ago.
At some point, I realised that the approach I was using could have a more general interest which, in my opinion, was worth a post.
In a few words, this post is about solving an algebra problem using a method originally designed to tackle statistical problems.

The Context
Recently @raj2018 submitted a question I'm going to resume this way:

Let S(phi ;  beta, f) a function of phi parameterized by beta and f.
Here is the graph of S(phi ;  0.449, 0.19)  @raj2018 provided

@raj2018 then asked how we can find other values (A, B)  of values for (beta, f) such that the graph of S(phi, A, B) has the same aspect of the graph above.
More precisely, let phi_0 the largest strictly negative value of phi such that  S(phi_0, A, B) = 0.
Then  S(phi, A, B) must be negative (strictly negative?) in the open interval (phi_0, 0), and must have exactly 3 extrema in this range.
I will said the point  (A, B) is admissible is S(phi, A, B) verifies thess conditions

The expression of S(phi, A, B) is that complex that it is likely impossible to find an (several?, all?) admissible point using analytic developments.

The approach

When I began thinking to this problem I early thought to find the entire domain of admissible points: was it something possible, at least with some reasonable accuracy? 

Quite rapidly I draw an analogy with an other type of problems whose solution is part of my job: the approximate construction of the probability density function (PDF) of multivariate random variables (obviously this implies that no analytical expression of this PDF is available). This is a very classical problem in Bayesian Statistics, for instance when we have to construt an approximation of a posterior PDF.

To stick with this example and put aside the rare situations where this PDF can be derived analytically, building a posterior PDF is largely based on specific numerical methods. 
The iconic one is known under the generic name MCMC  which stands for Markov Chain Monte Carlo.

Why am I speaking about MCMC or PDF or even random variables?
Let us consider some multivariate random variable R whose PDF as a constant on some bounded domain D and is equal to 0 elsewhere. R is then a uniform random variable with support supp(R) = D.
Assuming the domain Adm of admissible (beta, f) is bounded, we may  think of it as the support of some uniform random variable. Following this analogy we may expect to use some MCMC method to "build the PDF of the bivariate random variable (beta, f)", otherwise stated "to capture​​​​​​ the boundary of​ Adm".

The simplest MCMC method is the Metropolis-Hastings algorithm (MH).
In a few simple words MH builds a Markov chain this way:

1. Let us assume that the chain already contains elements e₁, ..., e_n.
   Let  f  some suitable "fitness" function (whose nature is of no importance right now).
2. A potential new element c is randomly picked in some neighborhood or e_n.
3. If the ratio f (c) / f (e_n) is larger than 1, we decide to put c into the chain (thus e_n+1 = c) otherwise we leave it to chance to decide whether or not c iis put into the chain.
   If chance decides the contrary,  then e_n is duclicated (thus e_n+1 = e_n).

MH is not the most efficient MCMC algorithm but it is efficient enough for what we want to achieve.
The main difficulty here is that there is no natural way to build the fitness function  f , mainly because the equivalent random variable I talked about is a purely abstract construction.

A preliminary observation is that if S(phi, beta, f) < 0 whatever phi in (phi_0, 0), then S has an odd number of extrema in (phi_0, 0). The simplest way to find these extrema is to search for the zeros of the derivative S' of S with respect to phi, while discardinq those where the second derivative can reveal "false" extrema where both S'' of S is null (I emphasize this last point because I didn't account for it in attached file).
The algorithm designed in this file probably misses a few points for not checking if S''=0, but it is important to keep in mind that we don't want a complete identification of  Adm but just the capture of its boundary.
Unless we are extremely unlucky there is only a very small chance that omitting to check if S''=0 will deeply modify this boundary.

How to define function f  ?
What we want is that  f (c) / f (e_n) represents the probability to decide wether c is an admissible point or not. In a Markov chain this  ratio represents how better or worse c is relatively to e_n, and this is essential for the chain to be a true Markov chain.
But as our aim is not to build a true Markov chain but simply a chain which looks like a Markov chain, we we can take some liberties and replace  f (c) / f (e_n) by some function  g(c) which quantifies the propability for c to be an admissible couple. So we want that  g(c) = 1 if  S(phi, c) has exactly M=3 negative extrema and  g(c) < 1 if M <> 3.
The previous algorihm transforms into:

1. Let us assume that the chain already contains elements e₁, ..., e_n.
   Let  g  a function which the propability that element is admissible
2. A potential new element c is randomly picked in some neighborhood or e_n.
3. If the ratio g(c) is larger than 1, we decide to put c into the chain (thus e_n+1 = c) otherwise we leave it to chance to decide whether or not c iis put into the chain.
   If chance decides the contrary,  then e_n is duclicated (thus e_n+1 = e_n).

This algorithm can also be seen as a kind of genetic algorithm.

A possible choice is  g(c)= exp(-|3-M|).
In the attached file I use instead the expression g(c) = (M + 1) / 4 fo several reasons:

• It is less sharp at M=3 and thus enables more often to put c into the chain, which increases its exploratory capabilities.
• The case M > 3, which no preliminary investigation was able to uncover, is by construction eliminated in the procedure Extrema which use an early stopping strategy (if as soon as more than M=3 negative extrema are found the procedure stops).

The algorithm I designed basically relies upon two stages:

1. The first one is aimed to construct a "long" Markov-like chain ("long" and not long because Markov chains are usually much longer than those I use).
   There are two goals here:
   1. Check if Adm is or not simply-connected or not (if it has holes or not).
   2. Find a first set of admissible points that can be used as starting points for subsequent chains.
       
2. Run several independent Markov-like chains from a reduced set of admissible points.
   The way this reduced set is constructed depends on the goal to achieve:
   1. One may think of adding points among those already known in order to assess the connectivity of Adm,
   2. or refinining the boundary of Adm.

These two concurent objectives are mixed in an ad hoc way depending on the observation of the results already in hand.

We point here an important feature of MCMC methods: behind their apparent algorithmic simplicity, it is common that high-quality results can only be obtained efficiently at the cost of problem-dependent tuning.

A last word to say that after several trials and failures I found it simpler to reparameterize the problems in terms of (phi_0, f) instead of (beta, f).

Codes and results

Choice g(c) = (M + 1) / 4 
The code : Extrema_and_MCMC.mw

To access the full results I got load this m file (do not bother its extension, Mapleprimes doesn't enable uploading m files) MCMC_20231209_160249.mw (save it and change it's extension in to m instead mw)

EDITED: choice  g(c)= exp(-|3-M|)
Here are the files contzining the code and the results:
Extrema_and_MCMC_g2.mw
MCMC_20231211_084053.mw

To ease the comparison of the two sets of results I used the same random seeds inn both codes.
Comparing the results got around the first admissible point is straightforward.
It's more complex for @raj2018's solution because the first step of the algorithim (drawing of a sibgle chain of length 1000) finds six times more admissible point with g(c)= exp(-|3-M|) than with g(c) = (M + 1) / 4.                                 

Plots of physical quantities has significantly improved with Maple 2022. The updated useunits option makes unit conversion errors in plots very unlikely. A lot of time is saved when creating plots of physical quantities where values and units must be correct.

One final source of user errors remains: The manual entry of incorrect units in labels.

Below is a way to avoid such errors by computing labels with units for three prevalent axis labeling schemes.

Other desireable labels are given as a suggestion for future plot label enhancements where plot commands could provide formating functionality.

The rendering on this website adds double brakets ⟦ ⟧ wherever units are used. You have to open the document to see how Maple renders.

Simple plot example: Solar irradiance in space (1) (2) This plot has inconsistent axis labeling: •  The vertical axis has units but no name •  The horizontal axis has a name and units but they are not easily distinguishable. Misinterpretation is possible. Due to the close spacing the label could be read as a product of the dimension "time squared" (the time t times hours h is of the dimension time squared). Or the reader confounds name and units. (The use of italic fonts for names and roman fonts for units might not be noticeable and is a convention that is not used everywhere.)   The above labeling should be improved for communication, documentation or publication purposes.   A quick attempt using strings and the options useuints and labels. Axes are now consistent and can be interpreted unambiguously. Formatting can still be improved.   Unfortunately, using the options useunits (for unit conversion) and labels this way introduces a new source of user error when labels are entered with the wrong units.   A way to address this and to ensure unit error-free plotting of expressions of physical quantities is the following:   Step1: Define two lists, one for the units to display and the other for the names to display (3) Step2: Compute labels from the lists This step avoids the labeling error: No manual entry of units in labels required. (4)   Dimensionless labels  Double brackets Parenthesis The axis values equal physical quantities divided by their units. The algebraic equation , for example, is physically speaking correct. Most functions of Maple can process dimensionless expression of the kind  if G is given with appropriate physical units. This way of using physical quantities is consistent with ISO 80000.   Used in Maple to enter units in 2D-Math input mode Can be confounded with functional notation. Units are therefore often written as a whole word (e.g. seconds instead of s).     The time to produce the above three plots was about 10 Minutes. The most part was spent to get the typesetting of the second and third plot correct.   What takes significant more time (more a question of hours when Typesetting is used for the first time) are   Labels with or  formatting.   This formatting might be preferred but is unfortunately again not free from user errors. (I would probably use it if there was a simple and safe way). (5)       Remarks •  For two reasons, I have not given an example with the often used square brackets [ ] because:          Maple uses square brackets already for lists and indexing purposes,     and ISO 80000 uses square brackets as an operator that extracts the unit from a physical quantity (e.g.       ). •  Adding a unit to each value at axes ticks would definitely be a nice labeling feature for simple units. •  Programmatically analyzing the units defined in list  above and converting them in a generic way to a typesetting structure is not possible with a few high-level commands.   •  For inline quotients like in , an additional backslash must be entered in 2D-Math: \/   (6)      This will not prevent evaluation to a normal quotient but the input can be used to create an atomic variable (select with mouse -> 2-D Math -> Atomic Variable) (7)      This makes labeling much easier as compared to typesetting commands (compare to the above statements). (8) In any case it is a good idea to read ?plot,typesetting before experimenting with typesetting.

Axes_with_unit_labels.mw

My personal preference is for dimensionless labels.

Note:

The solution to avoid labeling errors works only for Maple 2022 and higher.

Some plot commands do not support plotting with units, or they do not fully support it yet.