Download PDE_and_BC_Example_with_RootOf_eigenvalues.mw

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

The Maple help contains a nice example of a puzzle solver named alphametic,  see  ?Iterator,Permute.
The goal is to determine the distinct digits represented by letters which satisfy a given  equation. The provided solved example is:

        "16540*781 = 12904836 + 12904"
i.e.  m=1, a=6, p=5 etc.

It's worth studying the program (written as a module) because it includes a few useful techniques.
I would suggest the following exercise for a Maple (young) programmer.

1.  When solving the "classical" puzzle  "FORTY+TEN+TEN=SIXTY", instead of the correct answer "29786+850+850=31486",   you will see
"2978Y+850+850=3148Y".
Try to find what's going on and correct the code.

2. The solutions given by alphametic include numbers beginning with a zero.
Execute e.g. .
Modify the code to produce only standard numbers.

Here is a simple collection of Python scripts that allows Maple code to be embedded directly within a LaTeX document. The scripts will process the LaTeX document, calling Maple as required, making the Maple output directly accessible in the LaTeX document. The source, including extensive examples, can be found at https://github.com/leo-brewin/hybrid-latex. This  library also supports inclusion of active Mathematica, Matlab, Python and Cadabra code within LaTeX documents.

Here is a simple example (based on maple/example-01 in the GitHub site)

\documentclass[12pt]{mpllatex}

\begin{document}

\section*{Calculus}

\begin{maple}
   ans := diff(x*sin(x),x):                          # mpl (ans.501,ans)
   ans := eval(diff(x*sin(x),x),x=Pi/4):             # mpl (ans.502,ans)
   ans := int(2*sin(x)^2, x=a..b):                   # mpl (ans.503,ans)
   ans := int(2*exp(-x^2),x=0..infinity):            # mpl (ans.504,ans)
   ans := ''int(2*exp(-x^2),x=0..infinity)'':        # mpl (lhs.504,ans)
   ans := int(int(x^2 + y^2,  y=0..x),x=0..1):       # mpl (ans.505,ans)
   ans := ''int(int(x^2 + y^2,  y=0..x),x=0..1)'':   # mpl (lhs.505,ans)
\end{maple}

\begin{align*}
   &\mpl*{ans.501}\\
   &\mpl*{ans.502}\\
   &\mpl*{ans.503}\\
   \mpl{lhs.504}&=\Mpl{ans.504}\\
   \mpl{lhs.505}&=\Mpl{ans.505}
\end{align*}

\end{document}

and here is the corresponding output

example-01.pdf

You might recall this image being shared on social media some time ago.

Source: http://cvcl.mit.edu/hybrid_gallery/monroe_einstein.html

Look closely and you see Albert Einstein. However, if you move further away (or make the image smaller), you see Marilyn Monroe.

To create the image, the high spatial frequency data from an image of Albert Einstein was added to the low spatial frequency data from an image of Marilyn Monroe. This approach was pioneered by Oliva et al. (2006) and is influenced by the multiscale processing of human vision.

• When we view objects near us, we see fine detail (that is, higher spatial frequencies dominate).

• However, when we view objects at a distance, the broad outline has greater influence (that is, lower spatial frequencies dominate).

I thought I'd try to create a similar image in Maple (get the complete application here).

Here's an overview of the approach (as outlined in Oliva et al., 2006). I used different source images of Einstein and Monroe.

Let's start by loading some packages and defining a few procedures.

restart:
with(ImageTools):
with(SignalProcessing):

fft_shift := proc(M)
   local nRows, nCols, quad_1, quad_2, quad_3, quad_4, cRows, cCols;
   nRows, nCols := LinearAlgebra:-Dimensions(M):
   cRows, cCols := ceil(nRows/2), ceil(nCols/2):
   quad_1 := M[1..cRows,      1..cCols]:
   quad_2 := M[1..cRows,      cCols + 1..-1]:  
   quad_3 := M[cRows + 1..-1, cCols + 1..-1]:
   quad_4 := M[cRows + 1..-1, 1..cCols]:
   return <<quad_3, quad_2 |quad_4, quad_1>>:
end proc:

PowerSpectrum2D := proc(M)
   return sqrt~(abs~(M))
end proc:

gaussian_filter := (a, b, sigma) -> Matrix(2 * a, 2 * b, (i, j) -> evalf(exp(-((i - a)^2 + (j - b)^2) / (2 * sigma^2))), datatype = float[8]):

fft_shift() swaps quadrants of a 2D Fourier transform around so that the zero frequency components are in the center.

PowerSpectrum2D() returns the spectra of a 2D Fourier transform

gaussian_filter() will be used to apply a high or low-pass filter in the frequency domain (a and b are the number of rows and columns in the 2D Fourier transform, and sigma is the cut-off frequency.

Now let's import and display the original Einstein and Monroe images (both are the same size).

einstein_img := Read("einstein.png")[..,..,1]:
Embed(einstein_img)

marilyn_img  := Read("monroe.png")[..,..,1]:
Embed(monroe_img)

Let's convert both images to the spatial frequency domain (not many people know that SignalProcessing:-FFT can calculate the Fourier transform of matrices).

einstein_fourier := fft_shift(FFT(einstein_img)):
monroe_fourier   := fft_shift(FFT(monroe_img)):

Visualizing the power spectra of the unfiltered and filtered images isn't necessary, but helps illustrate what we're doing in the frequency domain.

First the spectra of the Einstein image. Lower frequency data is near the center, while higher frequency data is further away from the center.

Embed(Create(PowerSpectrum2D(einstein_fourier)))

Now the spectra of the Monroe image.

Embed(Create(PowerSpectrum2D(monroe_fourier)))

Now we need to filter the frequency content of both images.

First, define the cutoff frequencies for the high and low pass Gaussian filters.

sigma_einstein := 25:
sigma_monroe   := 10:

In the frequency domain, apply a high pass filter to the Einstein image, and a low pass filter to the Monroe image.

nRows, nCols := LinearAlgebra:-Dimension(einstein_img):

einstein_fourier_high_pass := einstein_fourier *~ (1 -~ gaussian_filter(nRows/2, nCols/2, sigma_einstein)):
monroe_fourier_low_pass    := monroe_fourier   *~ gaussian_filter(nRows/2, nCols/2, sigma_monroe):

Here's the spectra of the Einstein and Monroe images after the filtering (compare these to the pre-filtered spectra above).

Embed(Create(PowerSpectrum2D(einstein_fourier_high_pass)))

Embed(Create(PowerSpectrum2D(monroe_fourier_low_pass)))

Before combining both images in the Fourier domain, let's look the individual filtered images.

einstein_high_pass_img := Re~(InverseFFT(fft_shift(einstein_fourier_high_pass))):
monroe_low_pass_img    := Re~(InverseFFT(fft_shift(monroe_fourier_low_pass))):

We're left with sharp detail in the Einstein image.

Embed(FitIntensity(Create(einstein_high_pass_img)))

But the Monroe image is blurry, with only lower spatial frequency data.

Embed(FitIntensity(Create(monroe_low_pass_img)))

For the final image, we're simply going to add the Fourier transforms of both filtered images, and invert to the spatial domain.

hybrid_image := Create(Re~(InverseFFT(fft_shift(monroe_fourier_low_pass + einstein_fourier_high_pass)))):
Embed(hybrid_image)

So that's our final image, and has a similar property to the hybrid image at the top of this post.

• Move close to the computer monitor and you see Albert Einstein.
• Move to the other side of the room, and Marilyn Monroe swims into vision (if you're myopic, just take off your glasses and don't move back as much).

To simulate this, here, I've successively reduced the size of the hybrid image

And just because I can, here's a hybrid image of a cat and a dog, generated by the same worksheet.

I have recently visited the Queen's House at Greenwich  (see wiki),  an  important building in British architectural history (17th century).
I was impressed by the Great Hall Floor, whose central geometric decoration seems to be generated by a Maple program :-)

Here is my code for this. I hope you will like it.

restart;
with(plots): with(plottools):
n:=32: m:=3:#    n:=64: m:=7:

a[0], b[0] := exp(-Pi*I/n), exp(Pi*I/n):
c[0]:=b[0]+(a[0]-b[0])/sqrt(2)*exp(I*Pi/4):  
for k to m+1 do  
  c[k]:=a[k-1]+b[k-1]-c[k-1];
  b[k]:=c[k]*(1+exp(2*Pi*I/n))-b[k-1];
  a[k]:=conjugate(b[k])  od:
b[-1]:=c[0]*(1+exp(2*Pi*I/n))-b[0]:
a[-1]:=conjugate(b[-1]):
c[-1]:=a[-1]+b[-1]-c[0]:
seq( map[inplace](evalf@[Re,Im], w), w=[a,b,c] ):
Q:=polygonplot([seq([a[k],c[k],b[k],c[k+1]],k=0..m-1), [c[m],a[m],b[m]], [a[-1],b[-1],c[0]]]):
display(seq(rotate(Q, 2*k*Pi/n, [0,0]),k=0..n-1), disk([0,0],c[m][1]/3), axes=none, size=[800,800], color=black);

In this app you can visualize the location of the points in the different quadrants, also calculate the distance between two points. Finally the calculation of the coordinates of the midpoint. With these applications can be combined to study different cases between distance between two points and midpoint. Generated in Maple for students of secondary education and pre-calculation. In Spanish

Distance_between_two_points_and_midpoint.mw

Lenin Araujo Castillo

Ambassador of Maple

We produce only one rotation of the cube by the angle Pi/2, and then repeat it at the following points, changing the colors of the faces in turn. And so the illusion is created that the cube is rolling along a straight line without slipping.
(Just a picture without any sense.) cube_without_slipping.mw

Summary:Hough transform is a feature extraction algorithm widely used in digital image analysis. It identifies objects with specified shapes by letting all edge pixels in the image “vote” and extract candidates with the highest votes. The classical Hough transform is associated with the identification of straight lines, but with some modification, the algorithm can be used to detect objects of arbitrary shapes (usually circles, ellipses) in an image. Therefore, Hough transform has various applications in the field of computer vision. Samir and I implemented Hough transform to detect straight lines and circles in real world images and the results are below:

You can access the applications here for more details and try the algorithm on your own image!

Line Detection with Hough Transform

Circle Detection with Hough Transform

For people who want a bit more of the theory behind:

The voting procedure in Hough transform starts with an edge detector (we used the sobel edge detector in the applications). After locating all the edge pixels, the program will iterate through each of the pixels and record their “votes”. It is important to address that the voting procedure is completed in parameter space. For example, a straight line is usually in the form of y = a*x+b, the parameter space will be in terms of a and b. Hence y=2*x + 3 will be a point (2,3)  in the a-b space. However, notice that we cannot represent vertical lines with this system; it would take infinite memory to store a vertical line. Thus we use the polar coordinates system x*cos(θ) + y*sin(θ) - p = 0 instead, resulting in a θ -p space. Note that an image pixel (x₁,y₁) is represented by the intersection of all lines that satisfy x₁*cos(θ) + y₁*sin(θ) - p = 0, which in θ-p space is denoted as a sinusoidal curve. Hence the intersection of two different sinusoidal curves in θ-p space is essentially a unique, straight line that passes through two different image pixels. An intersection of 500 votes means that specific line passes through 500 edge pixels, implying that it is likely a solid contour line in the original image. In this way, by converting all edge pixels to sinusoidal curves in θ-p space, the algorithm is letting all pixels “vote”. In order to record the sinusoidal curves, we need an accumulator of θ * p dimensions. In the line detection application we used a matrix such that <= θ <= 180 and -sqrt(r² +c² ) <= p <= sqrt(r² +c² ), where r and c are the height (rows) and width (columns) of the image.

The voting procedure for circle detection is similar; the difference is in how we construct the accumulator to accommodate different parametrization for circles. Notice that a circle has three essential pieces of information, the x,y coordinates of centre and the radius. This suggests that our accumulator might have to be in three dimensions. In practice, circular Hough transform is usually performed with known or estimated radius, when the radius is not known, the most straightforward solution is indeed to test a range of possible radii with a 3-dimensional accumulator. However, this might be quite expensive, especially for images with high level of details.

Feel free to leave a comment if you have any question  :)  

Hello MaplePrime users,

Yesterday, I posted my MagicPuzzles package to the MapleCloud. It is a collection of tools I have written for manipulating, solving, and visualizing puzzles like Magic Squares and Magic Stars. Here's a sample solution for each:



For the Magic Square, the numbers on each horizontal and vertical line, along with the numbers on each of the two diagonals, add up to 65.

The inaugural version has separate sub-packages for:

• Magic Annulai (my own name)
• Magic Hexagons
• Magic Squares
• Magic Stars

Moreover, each sub-package contains these commands:

• Equations(), to return the linear equations for the variables based on the Magic Sum;
• Constraints(), to return the conditions that prevent redundant solutions found by reflections and rotations;
• VerifySolution(), to confirm if a list of numbers is a solution;
• EquivalentSolutions(), to determine solutions equivalent to a given solution;
• PrimarySolution(), which takes a solution and returns the associated primary solution;
• Reflection() and Rotation(), to reflect and rotate a solution; and
• Draw(), to provide a visualization (like the ones above).

There is also a command, MagicSolve(), which is used to find solutions, which take the form of permutations of [1,...n] for some positive integer n, to the equations. Essentially, it solves the linear equations, and cycles through all permutations for the free variables, and selects those that give "magic" solutions.

In future versions, I intend to add:

• Other specific classes of problems;
• More sample solutions; and
• Known algorithms for finding particular solutions.

To install the package, you can do so from here, or just execute the following within Maple 2017+:

PackageTools:-Install( 5755630338965504, 'overwrite' );

There are many examples in the help pages.

I think others will find this package interesting and useful, and I encourage you to check it out.

Allow me to introduce the Pauli Algebra package for Maple. This package implements the Clifford Algebra of Physical Space Cl₃ through the use of complex paravectors. The syntax of the package is similar to the notation popularized by the work of Dr. William E. Baylis. For more information, check out the Wikipedia entry on Paravectors and the APS workbook available on Dr. Baylis' website

To install the Pauli Algebra package while in a Maple session, type

PackageTools:-Install("https://maple.cloud/downloadDocument?id=5763496974221312&version=1")

Alternatively, you can find the Pauli Algebra package in the MapleCloud, and choose the install button on the far right.

Upon installation, several help files become available in Maple Help to assist you with the syntax.

As you will see in the Pauli.mpl source code of the package workbook, a number of Maple functions have been overloaded to handle the package's Paravector datatype. If you wish to overload additional Maple functions, let me know, and I'll include them in future updates.

I'll sign off by attaching a Maple worksheet/pdf containing some examples.

PauliAlgebraDemo.mw

PauliAlgebraDemo.pdf

Enjoy!

These files contain the kinematics and dynamics of the solid using a new technique (ALT + ENTER) to visualize the results online and thus save space in our Maple worksheet. Seen from a relative approach. For engineering students. In Spanish.

Equations_of_movement_for_particle.mw  (Intro)

Kinematics_and_relative_dynamics_of_a_solid.mw

Flat_kinetic_of_a_rigid_body.mw

Lenin Araujo Castillo

Ambassador of Maple

Recently I looked through an interesting book D. Wells "The Penquin book of Curious and Interesting Geometry" and came across this result, which I did not know about before: starting with a given quadrilateral , construct a square on each side. Van Aubel's theorem states that the two line segments between the centers of opposite squares are of equal lengths and are at right angles to one another. See the picture below:

It is interesting that this is true not only for a convex quadrilateral, but for arbitrary one and even self-intersecting one. This post is devoted to proving this result in Maple. The proof was very short and simple. Note that the coordinates of points are given not by lists, but by vectors. This is convenient when using  LinearAlgebra:-DotProduct  and  LinearAlgebra:-Norm  commands.

The code of the proof (the notation of the points on the picture coincide with their names in the code):

restart;
with(LinearAlgebra):
assign(seq(P||i=<x[i],y[i]>, i=1..4)):
P||5:=P||1:
assign(seq(Q||i=(P||i+P||(i+1))/2+<0,1; -1,0>.(P||(i+1)-P||i)/2, i=1..4)):
expand(DotProduct((Q||3-Q||1),(Q||4-Q||2), conjugate=false));
is(Norm(Q||3-Q||1, 2)=Norm(Q||4-Q||2, 2));

The output:

                                                      0
                                                    true

VA.mw

Hi MaplePrimes Users!

It’s your friendly, neighborhood tech support team; here to share some tips and tricks from issues we help users with on a daily basis.

A customer contacted us through a Help Page feedback form, asking how to add a row or column in a Matrix. The form came from the Row Operations help page, but the wording of the message suggested that the customer actually wanted to insert a new row or column altogether. Such manipulations can often be accomplished by a command in the ArrayTools package, but the only Insert command currently available is the one for Vectors and 1-D Arrays. Using the Concatenate command from that package, and various commands from the LinearAlgebra package (including the SubMatrix command), we were able to write two custom procedures to perform these manipulations:

InsertRow := proc (A::rtable, n::integer, v::Vector[row])
    local R, C, top, bottom;
    uses LinearAlgebra;
    R := RowDimension(A); C := ColumnDimension(A);
    top := SubMatrix(A, [1 .. n-1], [1 .. C]);
    bottom := SubMatrix(A, [n .. R], [1 .. C]);
    return ArrayTools:-Concatenate(1, top, v, bottom);
end proc:

InsertColumn := proc (A::rtable, n::integer, v::Vector[column])
    local R, C, left, right;
    uses LinearAlgebra;
    R := RowDimension(A); C := ColumnDimension(A);
    left := SubMatrix(A, [1 .. R], [1 .. n-1]);
    right := SubMatrix(A, [1 .. R], [n .. C]);
    return ArrayTools:-Concatenate(2, left, v, right)
end proc:

# test cases:

M := Matrix([[1, 2, 3], [4, 5, 6], [7, 8, 9]]):
v := Vector[row]([2, 2, 2]):
v2 := Vector[column]([2, 2, 2]):
seq(InsertRow(M, i, v), i = 1 .. 4);
seq(InsertColumn(M, i, v2), i = 1 .. 4);

We then reworked this problem using some handy indexing and construction notation that allows our previous procedures to save on the variable constructions and syntax:

InsertRow := proc( A :: rtable, V :: Vector[row], r :: posint )
    return < A[1..r-1,..]; V; A[r..-1,..] >:
end proc:

InsertColumn := proc( A :: rtable, V :: Vector[column], c :: posint )
    return < A[..,1..c-1] | V | A[..,c..-1] >:
end proc:

M := Matrix(3, 3, [seq(i, i = 1 .. 9)]);
A := convert(M, Array);
U := Vector[row]( [ a, b, c ] );
V := convert( U, 'Vector[column]' );
seq(InsertRow( A, U, i ), i=1..4);
seq(InsertColumn( A, V, i ), i=1..4);
seq(InsertRow( M, U, i ), i=1..4);
seq(InsertColumn( M, V, i ), i=1..4);

In order to explore the use of signal processing package in image processing, @Samir Khan and I created this application that performs discrete Haar wavelet transform on images to achieve both lossy (irreversible) and lossless (reversible) compression.

Haar wavelet compression modifies the image matrix to increase the number of zero entries in the matrix, which results in a sparse matrix that can be stored more efficiently, thus reducing the file size. A Haar wavelet transform groups adjacent items in the matrix, stores the average and difference of the two numbers. Notice that this computation is reversible since knowing the values of a, b and given that x1-x2 = a; (x1+x2)/2 = b, we can solve for x1 and x2. Haar wavelet compression is taking advantage of the property that neighboring pixels in an image usually share very similar value; hence recursively applying Haar wavelet transform to the rows and columns of an image matrix significantly increases the number of zero entries. In the application we achieved a compression ratio of 1.296 (number of non-zero entries in original: number of non-zero entries in processed matrix).

The fact that Haar wavelet transform is reversible means that we can use it to perform lossless image compression: the decompressed image will be exactly the same as the image before compression. Transmission and temporary storage of the data would be made more efficient without any loss of details in the image.

But what if the file size is still too big or we simply don’t need that many details in the image? That is when lossy compression comes into use. By omitting some details/fidelity, lossy compression is able to achieve notably smaller file size. In this application, we apply a filter to the transformed image matrix, setting entries that are close to zeros to actual zeros (i.e. pick a positive number ϵ such that all x < ϵ are changed to 0 in the matrix). The value of ϵ directly impacts the quality of the decompressed image so should be chosen carefully in practice. In this application, we chose ϵ = 0.01, which results in a compression ratio of 19.38, but instead produces a very blurry image after reversing the compression.

(left: Original image, right: lossy compression with ϵ = 0.01)

The application can be accessed here for more details.

Download What_is_New_after_Maple_2018.mw

Katherina von Bülow