Education

Teaching and learning about math, Maple and MapleSim

Here at Maplesoft, we like to foster innovation in technological development. Whether that is finding solutions to global warming, making medical discoveries that save millions, or introducing society to very advanced functional robots, Maplesoft is happy to contribute, support and encourage innovative people and organizations researching these complex topics. This year, we are delighted to have sponsored two contests in the robotics field that provide opportunities to think big and make an impact: Create the Future Design Contest and the International Space Apps Challenge. 

Create the Future Design Contest

Established in 2002, and organized by TechBriefs, the goal of the Create the Future Design Contest is to help engineers bring their product design ideas to life. The overall ‘mission of the contest is to benefit humanity, the environment, and the economy.’ This year, there were a record 1,159 new product ideas submitted by students, engineers, and entrepreneurs from all over the globe. In the machinery/automation/robotics category, which Maplesoft sponsored, the project with the top votes was designed by two engineers who chose to name their innovation CAP Exoskeleton, a type of assistive robotic machine designed to aid the user in walking, squatting, and carrying heavy loads over considerable distances. It can either be used to enhance physical endurance for military purposes or to help the physically impaired perform daily tasks. A contest like Create the Future is a perfect opportunity, for engineers in particular, to learn, explore, and create. 

The CAP Exoskeleton - ©2015 Create the Future Design Contest

 

International Space Apps Challenge

The exploration of space has always been unique in its search for knowledge. The International Space Apps Challenge, a NASA incubator innovation program, is an ‘international mass collaboration focused on space exploration that takes place over 48-hours in cities around the world’. It is a unique global competition where people rally together to find solutions to real world problems, bringing humanity closer to understanding the Earth, the universe, the human race, and robotics. These goals, the organizers believe, can be reached much faster if we combine the power of the seven billion or so brains that occupy the planet, not forgetting the six that are currently orbiting above us aboard the International Space Station. The competition is open to people of all ages and in all fields, including engineers, technologists, scientists, designers, artists, educators, students, entrepreneurs, and so on. With an astounding 13,846 participants from all over the world, several highly innovative solutions were presented. 

Maplesoft sponsored the University of York location in the UK where the winning team of five modeled an app called CropOp, a communication tool that connects the government to local farmers with the goal of providing instantaneous, crucial information regarding pest breakout warnings, extreme weather, and other important updates. This UK-based team believes the quality and quantity of food produced will be improved, especially benefiting the undernourished communities in Africa. Maplesoft supports the Space Apps Challenge because it proves that collaboration makes for bigger and better discoveries that can save millions of people.

 

Donating Maplesoft software for contestants to use is part of the sponsorship. The real delight is to wait and see what innovative concepts they come up with. When we sponsor contests like these, we find it benefits our software as much as it does the participants. Plus, if the contestants can provide solutions to real world issues, well, that benefits everyone! 

ABSTRACT. In this paper we demonstrate how the simulation of dynamic systems engineering has been implemented with graphics software algorithms using maple and MapleSim. Today, many of our researchers the computational modeling performed by inserting a piece of code from static work; with these packages we have implemented through the automation components of kinematics and dynamics of solids simple to complex.

It is very important to note that once developed equations study; recently we can move to the simulation; to thereby start the physical construction of the system. We will use mathematical and computational methods using the embedded buttons which lie in the dynamics leaves and viewing platform cloud of Maplesoft and power MapleNet for online evaluation of specialists in the area. Finally they will see some work done; which integrate various mechanical and computational concepts implemented for companies in real time and pattern of credibility.

 

Selasi_2015.pdf

(in spanish)

 

Lenin Araujo Castillo

 

 

I have two linear algebra texts [1, 2]  with examples of the process of constructing the transition matrix Q that brings a matrix A to its Jordan form J. In each, the authors make what seems to be arbitrary selections of basis vectors via processes that do not seem algorithmic. So recently, while looking at some other calculations in linear algebra, I decided to revisit these calculations in as orderly a way as possible.

 

First, I needed a matrix A with a prescribed Jordan form. Actually, I started with a Jordan form, and then constructed A via a similarity transform on J. To avoid introducing fractions, I sought transition matrices P with determinant 1.

 

Let's begin with J, obtained with Maple's JordanBlockMatrix command.

 

• 

Tools_Load Package: Linear Algebra

Loading LinearAlgebra

J := JordanBlockMatrix([[2, 3], [2, 2], [2, 1]])

Matrix([[2, 1, 0, 0, 0, 0], [0, 2, 1, 0, 0, 0], [0, 0, 2, 0, 0, 0], [0, 0, 0, 2, 1, 0], [0, 0, 0, 0, 2, 0], [0, 0, 0, 0, 0, 2]])

 

``

The eigenvalue lambda = 2 has algebraic multiplicity 6. There are sub-blocks of size 3×3, 2×2, and 1×1. Consequently, there will be three eigenvectors, supporting chains of generalized eigenvectors having total lengths 3, 2, and 1. Before delving further into structural theory, we next find a transition matrix P with which to fabricate A = P*J*(1/P).

 

The following code generates random 6×6 matrices of determinant 1, and with integer entries in the interval [-2, 2]. For each, the matrix A = P*J*(1/P) is computed. From these candidates, one A is then chosen.

 

L := NULL:

 

 

After several such trials, the matrix A was chosen as

 

A := Matrix(6, 6, {(1, 1) = -8, (1, 2) = -8, (1, 3) = 4, (1, 4) = -8, (1, 5) = -1, (1, 6) = 5, (2, 1) = -1, (2, 2) = 3, (2, 3) = 1, (2, 4) = -2, (2, 5) = 2, (2, 6) = -1, (3, 1) = -13, (3, 2) = -9, (3, 3) = 8, (3, 4) = -11, (3, 5) = 1, (3, 6) = 5, (4, 1) = 3, (4, 2) = 3, (4, 3) = -1, (4, 4) = 4, (4, 5) = 1, (4, 6) = -2, (5, 1) = 7, (5, 2) = 5, (5, 3) = -3, (5, 4) = 6, (5, 5) = 2, (5, 6) = -3, (6, 1) = -6, (6, 2) = -2, (6, 3) = 3, (6, 4) = -7, (6, 5) = 2, (6, 6) = 3})

 

 

for which the characteristic and minimal polynomials are

 

factor(CharacteristicPolynomial(A, lambda))

(lambda-2)^6

factor(MinimalPolynomial(A, lambda))

(lambda-2)^3

 

 

So, if we had started with just A, we'd now know that the algebraic multiplicity of its one eigenvalue lambda = 2 is 6, and there is at least one 3×3 sub-block in the Jordan form. We would not know if the other sub-blocks were all 1×1, or a 1×1 and a 2×2, or another 3×3. Here is where some additional theory must be invoked.

``

The null spaces M[k] of the matrices (A-2*I)^k are nested: `⊂`(`⊂`(M[1], M[2]), M[3]) .. (), as depicted in Figure 1, where the vectors a[k], k = 1, () .. (), 6, are basis vectors.

 

Figure 1   The nesting of the null spaces M[k] 

 

 

The vectors a[1], a[2], a[3] are eigenvectors, and form a basis for the eigenspace M[1]. The vectors a[k], k = 1, () .. (), 5, form a basis for the subspace M[2], and the vectors a[k], k = 1, () .. (), 6, for a basis for the space M[3], but the vectors a[4], a[5], a[6] are not yet the generalized eigenvectors. The vector a[6] must be replaced with a vector b[6] that lies in M[3] but is not in M[2]. Once such a vector is found, then a[4] can be replaced with the generalized eigenvector `≡`(b[4], (A-2*I)^2)*b[6], and a[1] can be replaced with `≡`(b[1], A-2*I)*b[4]. The vectors b[1], b[4], b[6] are then said to form a chain, with b[1] being the eigenvector, and b[4] and b[6] being the generalized eigenvectors.

 

If we could carry out these steps, we'd be in the state depicted in Figure 2.

 

Figure 2   The null spaces M[k] with the longest chain determined

 

 

Next, basis vector a[5] is to be replaced with b[5], a vector in M[2] but not in M[1], and linearly independent of b[4]. If such a b[5] is found, then a[2] is replaced with the generalized eigenvector `≡`(b[2], A-2*I)*b[5]. The vectors b[2] and b[5] would form a second chain, with b[2] as the eigenvector, and b[5] as the generalized eigenvector.

``

Define the matrix C = A-2*I by the Maple calculation

 

C := A-2

Matrix([[-10, -8, 4, -8, -1, 5], [-1, 1, 1, -2, 2, -1], [-13, -9, 6, -11, 1, 5], [3, 3, -1, 2, 1, -2], [7, 5, -3, 6, 0, -3], [-6, -2, 3, -7, 2, 1]])

 

``

and note

 

N := convert(NullSpace(C), list)

[Vector(6, {(1) = 1/2, (2) = 1/2, (3) = 1, (4) = 0, (5) = 0, (6) = 1}), Vector(6, {(1) = -1/2, (2) = -1/2, (3) = -2, (4) = 0, (5) = 1, (6) = 0}), Vector(6, {(1) = -2, (2) = 1, (3) = -1, (4) = 1, (5) = 0, (6) = 0})]

NN := convert(LinearAlgebra:-NullSpace(C^2), list)

[Vector(6, {(1) = 2/5, (2) = 0, (3) = 0, (4) = 0, (5) = 0, (6) = 1}), Vector(6, {(1) = 0, (2) = 0, (3) = 0, (4) = 0, (5) = 1, (6) = 0}), Vector(6, {(1) = -1, (2) = 0, (3) = 0, (4) = 1, (5) = 0, (6) = 0}), Vector(6, {(1) = 2/5, (2) = 0, (3) = 1, (4) = 0, (5) = 0, (6) = 0}), Vector(6, {(1) = -3/5, (2) = 1, (3) = 0, (4) = 0, (5) = 0, (6) = 0})]

 

``

The dimension of M[1] is 3, and of M[2], 5. However, the basis vectors Maple has chosen for M[2] do not include the exact basis vectors chosen for M[1].

 

We now come to the crucial step, finding b[6], a vector in M[3] that is not in M[2] (and consequently, not in M[1] either). The examples in [1, 2] are simple enough that the authors can "guess" at the vector to be taken as b[6]. What we will do is take an arbitrary vector in M[3] and project it onto the 5-dimensional subspace M[2], and take the orthogonal complement as b[6].

``

A general vector in M[3] is

 

Z := `<,>`(u || (1 .. 6))

Vector[column]([[u1], [u2], [u3], [u4], [u5], [u6]])

 

``

A matrix that projects onto M[2] is

 

P := ProjectionMatrix(NN)

Matrix([[42/67, -15/67, 10/67, -25/67, 0, 10/67], [-15/67, 58/67, 6/67, -15/67, 0, 6/67], [10/67, 6/67, 63/67, 10/67, 0, -4/67], [-25/67, -15/67, 10/67, 42/67, 0, 10/67], [0, 0, 0, 0, 1, 0], [10/67, 6/67, -4/67, 10/67, 0, 63/67]])

 

``

The orthogonal complement of the projection of Z onto M[2] is then -P*Z+Z. This vector can be simplified by choosing the parameters in Z appropriately. The result is taken as b[6].

 

b[6] := 67*(eval(Z-Typesetting:-delayDotProduct(P, Z), Equate(Z, UnitVector(1, 6))))*(1/5)

Vector[column]([[5], [3], [-2], [5], [0], [-2]])

NULL

 

``

The other two members of this chain are then

 

b[4] := Typesetting:-delayDotProduct(C, b[6])

Vector[column]([[-132], [-12], [-169], [40], [92], [-79]])

b[1] := Typesetting:-delayDotProduct(C, b[4])

Vector[column]([[-67], [134], [67], [67], [0], [134]])

 

``

A general vector in M[2] is a linear combination of the five vectors that span the null space of C^2, namely, the vectors in the list NN. We obtain this vector as

 

ZZ := add(u || k*NN[k], k = 1 .. 5)

Vector[column]([[(2/5)*u1-u3+(2/5)*u4-(3/5)*u5], [u5], [u4], [u3], [u2], [u1]])

 

``

A vector in M[2] that is not in M[1] is the orthogonal complement of the projection of ZZ onto the space spanned by the eigenvectors spanning M[1] and the vector b[4]. This projection matrix is

 

PP := LinearAlgebra:-ProjectionMatrix(convert(`union`(LinearAlgebra:-NullSpace(C), {b[4]}), list))

Matrix([[69/112, -33/112, 19/112, -17/56, 0, 19/112], [-33/112, 45/112, 25/112, 13/56, 0, 25/112], [19/112, 25/112, 101/112, 1/56, 0, -11/112], [-17/56, 13/56, 1/56, 5/28, 0, 1/56], [0, 0, 0, 0, 1, 0], [19/112, 25/112, -11/112, 1/56, 0, 101/112]])

 

``

The orthogonal complement of ZZ, taken as b[5], is then

 

b[5] := 560*(eval(ZZ-Typesetting:-delayDotProduct(PP, ZZ), Equate(`<,>`(u || (1 .. 5)), LinearAlgebra:-UnitVector(4, 5))))

Vector[column]([[-9], [-59], [17], [58], [0], [17]])

 

``

Replace the vector a[2] with b[2], obtained as

 

b[2] := Typesetting:-delayDotProduct(C, b[5])

Vector[column]([[251], [-166], [197], [-139], [-112], [-166]])

 

 

The columns of the transition matrix Q can be taken as the vectors b[1], b[4], b[6], b[2], b[5], and the eigenvector a[3]. Hence, Q is the matrix

 

Q := `<|>`(b[1], b[4], b[6], b[2], b[5], N[3])

Matrix([[-67, -132, 5, 251, -9, -2], [134, -12, 3, -166, -59, 1], [67, -169, -2, 197, 17, -1], [67, 40, 5, -139, 58, 1], [0, 92, 0, -112, 0, 0], [134, -79, -2, -166, 17, 0]])

 

``

Proof that this matrix Q indeed sends A to its Jordan form consists in the calculation

 

1/Q.A.Q = Matrix([[2, 1, 0, 0, 0, 0], [0, 2, 1, 0, 0, 0], [0, 0, 2, 0, 0, 0], [0, 0, 0, 2, 1, 0], [0, 0, 0, 0, 2, 0], [0, 0, 0, 0, 0, 2]])``

 

NULL

The bases for M[k], k = 1, 2, 3, are not unique. The columns of the matrix Q provide one set of basis vectors, but the columns of the transition matrix generated by Maple, shown below, provide another.

 

JordanForm(A, output = 'Q')

Matrix([[-5, -43/5, -9/5, 7/5, -14/5, -3/5], [10, -4/5, -6/25, 1/5, -6/25, -3/25], [5, -52/5, -78/25, 13/5, -78/25, -39/25], [5, 13/5, 38/25, -2/5, 38/25, 4/25], [0, 6, 42/25, -1, 42/25, 21/25], [10, -29/5, -11/25, 1/5, -11/25, 7/25]])

 

``

I've therefore added to my to-do list the investigation into Maple's algorithm for determining an appropriate set of basis vectors that will support the Jordan form of a matrix.

 

References

 

NULL

[1] Linear Algebra and Matrix Theory, Evar Nering, John Wiley and Sons, Inc., 1963

[2] Matrix Methods: An Introduction, Richard Bronson, Academic Press, 1969

 

NULL

``

Download JordanForm_blog.mw

This October 21st, Maplesoft will be hosting a full-production, live streaming webinar featuring Dr. Robert Lopez, Emeritus Professor and Maple Fellow. You might have caught Dr. Lopez's Clickable Calculus webinar series before, but this webinar is your chance to meet the man behind the voice and watch him use Clickable Math techniques live!

In this webinar, Dr. Lopez will present examples of what "resequencing concepts and skills" looks like when implemented with Maple's point-and-click syntax-free paradigm. He will demonstrate how Maple can not only be used to elucidate the concept, but also, how it can be used to illustrate and implement the manipulations that ultimately the student must master.

Click here for more information and registration.

In this paper we will demonstrate the many differences of implementation in the modeling of mechanical systems using embedded components through Maplesoft. The mechanical systems are used for different tasks and therefore have different structure in its design; as to the nature of the used functional elements placed on them, they vary greatly. This diversity is reflected in approaches and practices in modeling.

The following cases focus on mechanical components of the units manufacturing and processing machines. We can generate graphs for analysis using different dynamic pair ametros; all in real-time considerations in its manufacturing costs from the equations of conservation of energy.
Therefore modeling with Maplesoft ensures the smooth optimum performance in mechanical systems, highlighting the sustainability criteria for other areas of engineering.

 

XXXIII_Coloquio_SMP_2015.pdf

XXXIII_Coloquio_UNASAM_2015.mw

(in spanish)

L.AraujoC.

 

 


One of the interesting things about the Physics package is that it was designed from scratch to extend the domain of operations of the Maple system from commutative variables to one that includes commutative, anticommutative and nonocommutative variables, as well as abstract vectors and related (nabla) differential operators. In this line we have, among others, the following Physics commands working with this extended domain: `*` , `.` , `^` , diff , Expand , Normal , Simplify , Gtaylor , and Coefficients .

 

More recently, Pascal Szriftgiser (from Laboratoire PhLAM, Université Lille 1, France), suggested a similar approach to factorize expressions involving noncommutative variables. This is a pretty complicated problem though. Pascal's suggestion, however, spinned around an idea beautiful for its simplicity, similar to what is done in the experimental Physics command, PerformOnAnticommutativeSystem , that is, to transform the problem into one that can be treated with the command that works only with commutative variables and from there extract the result for noncommutative ones.The approach has limitations but it is surprising how far one can go using imaginative algebraic manipulations to extend these commands that otherwise only work with commutative variables.

 

In brief, we now have a new command, Physics:-Factor, with already powerful performance for factorizing algebraic expressions that involve commutative, noncommutative and anticommutative variables, making Maple's mathematical capabilities more advanced in very interesting directions. This command is in fact useful not just in advanced theoretical physics, but for instance also when working with noncommutative symbols representing abstract matrices (that can have dependency, and so they can be differentiated before saying anything about their components, multiplied, and be present int  expressions that in turn can be expanded, simplified and now also factorized), and also useful with expressions that include differential operators, now that within Physics you can compute with the the covariant and noncovariant derivatives D_  and d_ algebraically. For instance, how about solving differential equations using Physics:-Factor (reducing their order by means of factoring the involved differential operators) ? :)

 

What follows are some basic algebraic examples illustrating the novelty, and as usual to reproduce the results in this worksheet you need to update your Physics library with the one available in the Maplesoft R&D Physics webpage.

 

Physics:-Version()[2]

`2015, September 25, 7:48 hours`

(1)

with(Physics); -1; Setup(quantumoperators = {a, b, c, d, e}, mathematicalnotation = true)

[mathematicalnotation = true, quantumoperators = {a, b, c, d, e}]

(2)

First example, because of using mathematical notation, noncommutative variables are displayed in different color (olive)

Physics:-`*`(Physics:-`^`(alpha, 2), Physics:-`^`(a, 2))+Physics:-`*`(Physics:-`*`(Physics:-`*`(alpha, sqrt(2)), a), b)+Physics:-`*`(Physics:-`*`(Physics:-`*`(Physics:-`*`(4, sqrt(2)), lambda), Physics:-`^`(b, 2)), c)+Physics:-`*`(Physics:-`*`(Physics:-`*`(Physics:-`*`(Physics:-`*`(4, lambda), alpha), b), c), a)+Physics:-`*`(Physics:-`*`(Physics:-`*`(Physics:-`*`(Physics:-`*`(4, lambda), sqrt(2)), b), c), b)+Physics:-`*`(Physics:-`*`(16, Physics:-`^`(lambda, 2)), Physics:-`^`(Physics:-`*`(b, c), 2))+Physics:-`*`(Physics:-`*`(Physics:-`*`(Physics:-`*`(Physics:-`*`(4, alpha), lambda), a), b), c)+Physics:-`*`(Physics:-`*`(Physics:-`*`(sqrt(2), alpha), b), a)+Physics:-`*`(2, Physics:-`^`(b, 2))

alpha^2*Physics:-`^`(a, 2)+alpha*2^(1/2)*Physics:-`*`(a, b)+4*2^(1/2)*lambda*Physics:-`*`(Physics:-`^`(b, 2), c)+4*lambda*alpha*Physics:-`*`(b, c, a)+4*2^(1/2)*lambda*Physics:-`*`(b, c, b)+16*lambda^2*Physics:-`^`(Physics:-`*`(b, c), 2)+4*lambda*alpha*Physics:-`*`(a, b, c)+alpha*2^(1/2)*Physics:-`*`(b, a)+2*Physics:-`^`(b, 2)

(3)

Physics:-Factor(alpha^2*Physics:-`^`(a, 2)+alpha*2^(1/2)*Physics:-`*`(a, b)+4*2^(1/2)*lambda*Physics:-`*`(Physics:-`^`(b, 2), c)+4*lambda*alpha*Physics:-`*`(b, c, a)+4*2^(1/2)*lambda*Physics:-`*`(b, c, b)+16*lambda^2*Physics:-`^`(Physics:-`*`(b, c), 2)+4*lambda*alpha*Physics:-`*`(a, b, c)+alpha*2^(1/2)*Physics:-`*`(b, a)+2*Physics:-`^`(b, 2))

Physics:-`^`(4*lambda*Physics:-`*`(b, c)+a*alpha+2^(1/2)*b, 2)

(4)

A more involved example from a physics problem, illustrating that the factorization is also happening within function's arguments, as well as that we can also correctly expand mathematical expressions involving noncommutative variables

PDEtools:-declare((a, b, c, g)(x, y)):

a(x, y)*`will now be displayed as`*a

 

b(x, y)*`will now be displayed as`*b

 

c(x, y)*`will now be displayed as`*c

 

g(x, y)*`will now be displayed as`*g

(5)

Physics:-Intc(Physics:-`^`(Physics:-`*`(Physics:-`*`(Physics:-`*`(4, Physics:-Dagger(b(x, y))), c(x, y)), lambda)+Physics:-`*`(Physics:-`*`(Physics:-`*`(alpha, f(t)), a(x, y)), Physics:-Dagger(a(x, y)))+Physics:-`*`(Physics:-`*`(sqrt(2), g(x, y)), b(x, y)), 2), x, y)

Int(Int(Physics:-`^`(4*lambda*Physics:-`*`(Physics:-Dagger(b(x, y)), c(x, y))+alpha*f(t)*Physics:-`*`(a(x, y), Physics:-Dagger(a(x, y)))+2^(1/2)*g(x, y)*b(x, y), 2), x = -infinity .. infinity), y = -infinity .. infinity)

(6)

So first expand ...

expand(Int(Int(Physics:-`^`(4*lambda*Physics:-`*`(Physics:-Dagger(b(x, y)), c(x, y))+alpha*f(t)*Physics:-`*`(a(x, y), Physics:-Dagger(a(x, y)))+2^(1/2)*g(x, y)*b(x, y), 2), x = -infinity .. infinity), y = -infinity .. infinity))

Int(Int(16*lambda^2*Physics:-`*`(Physics:-Dagger(b(x, y)), c(x, y), Physics:-Dagger(b(x, y)), c(x, y))+4*lambda*alpha*f(t)*Physics:-`*`(Physics:-Dagger(b(x, y)), c(x, y), a(x, y), Physics:-Dagger(a(x, y)))+4*lambda*2^(1/2)*g(x, y)*Physics:-`*`(Physics:-Dagger(b(x, y)), c(x, y), b(x, y))+4*alpha*f(t)*lambda*Physics:-`*`(a(x, y), Physics:-Dagger(a(x, y)), Physics:-Dagger(b(x, y)), c(x, y))+alpha^2*f(t)^2*Physics:-`*`(a(x, y), Physics:-Dagger(a(x, y)), a(x, y), Physics:-Dagger(a(x, y)))+alpha*f(t)*2^(1/2)*g(x, y)*Physics:-`*`(a(x, y), Physics:-Dagger(a(x, y)), b(x, y))+4*2^(1/2)*g(x, y)*lambda*Physics:-`*`(b(x, y), Physics:-Dagger(b(x, y)), c(x, y))+2^(1/2)*g(x, y)*alpha*f(t)*Physics:-`*`(b(x, y), a(x, y), Physics:-Dagger(a(x, y)))+2*g(x, y)^2*Physics:-`^`(b(x, y), 2), x = -infinity .. infinity), y = -infinity .. infinity)

(7)

Now retrieve the original expression by recursing over the arguments and so factoring the integrand

Physics:-Factor(Int(Int(16*lambda^2*Physics:-`*`(Physics:-Dagger(b(x, y)), c(x, y), Physics:-Dagger(b(x, y)), c(x, y))+4*lambda*alpha*f(t)*Physics:-`*`(Physics:-Dagger(b(x, y)), c(x, y), a(x, y), Physics:-Dagger(a(x, y)))+4*lambda*2^(1/2)*g(x, y)*Physics:-`*`(Physics:-Dagger(b(x, y)), c(x, y), b(x, y))+4*alpha*f(t)*lambda*Physics:-`*`(a(x, y), Physics:-Dagger(a(x, y)), Physics:-Dagger(b(x, y)), c(x, y))+alpha^2*f(t)^2*Physics:-`*`(a(x, y), Physics:-Dagger(a(x, y)), a(x, y), Physics:-Dagger(a(x, y)))+alpha*f(t)*2^(1/2)*g(x, y)*Physics:-`*`(a(x, y), Physics:-Dagger(a(x, y)), b(x, y))+4*2^(1/2)*g(x, y)*lambda*Physics:-`*`(b(x, y), Physics:-Dagger(b(x, y)), c(x, y))+2^(1/2)*g(x, y)*alpha*f(t)*Physics:-`*`(b(x, y), a(x, y), Physics:-Dagger(a(x, y)))+2*g(x, y)^2*Physics:-`^`(b(x, y), 2), x = -infinity .. infinity), y = -infinity .. infinity))

Int(Int(Physics:-`^`(4*lambda*Physics:-`*`(Physics:-Dagger(b(x, y)), c(x, y))+alpha*f(t)*Physics:-`*`(a(x, y), Physics:-Dagger(a(x, y)))+2^(1/2)*g(x, y)*b(x, y), 2), x = -infinity .. infinity), y = -infinity .. infinity)

(8)

This following one looks simpler but it is actually more complicated:

Physics:-`*`(Physics:-Commutator(a, b), c)

Physics:-`*`(Physics:-Commutator(a, b), c)

(9)

expand(Physics:-`*`(Physics:-Commutator(a, b), c))

Physics:-`*`(a, b, c)-Physics:-`*`(b, a, c)

(10)

The complication consists of the fact that the standard factor  command, that assumes products are commutative, can never deal with factors like Physics:-Commutator(a, b) = a*b-a*b because if products were commutative these factors are equal to 0. Of course we not just us factor but include a number of algebraic manipulations before using it, so that the approach handles these cases nicely anyway

Physics:-Factor(Physics:-`*`(a, b, c)-Physics:-`*`(b, a, c))

Physics:-`*`(Physics:-`*`(a, b)-Physics:-`*`(b, a), c)

(11)

This other one is more complicated:

Physics:-`*`(Physics:-`*`(a, b)-Physics:-`*`(b, a), a+Physics:-`*`(beta, b)+Physics:-`^`(c, 2))

Physics:-`*`(Physics:-`*`(a, b)-Physics:-`*`(b, a), a+beta*b+Physics:-`^`(c, 2))

(12)

When you expand,

expand(Physics:-`*`(Physics:-`*`(a, b)-Physics:-`*`(b, a), a+beta*b+Physics:-`^`(c, 2)))

Physics:-`*`(a, b, a)+beta*Physics:-`*`(a, Physics:-`^`(b, 2))+Physics:-`*`(a, b, Physics:-`^`(c, 2))-Physics:-`*`(b, Physics:-`^`(a, 2))-beta*Physics:-`*`(b, a, b)-Physics:-`*`(b, a, Physics:-`^`(c, 2))

(13)

you see that there are various terms involving the same noncommutative operands, just multiplied in different order. Generally speaking the limitation (n this moment) of the approach is: "there cannot be more than 2 terms in the expanded form containing the same operands" . For instance in (13) the 1st and 4th terms have the same operands, that are actually also present in the 5th term but there you also have beta and for that reason (involving some additional manipulations) it can be handled:

Physics:-Factor(Physics:-`*`(a, b, a)+beta*Physics:-`*`(a, Physics:-`^`(b, 2))+Physics:-`*`(a, b, Physics:-`^`(c, 2))-Physics:-`*`(b, Physics:-`^`(a, 2))-beta*Physics:-`*`(b, a, b)-Physics:-`*`(b, a, Physics:-`^`(c, 2)))

Physics:-`*`(Physics:-`*`(a, b)-Physics:-`*`(b, a), a+beta*b+Physics:-`^`(c, 2))

(14)

Recalling, in all these examples, the task is actually accomplished by the standard factor  command, and the manipulations consist of ingeniously rewriting the given problem as one that involves only commutative variables, and from extract the correct result for non commutative variables.

 

To conclude, here is an example where the approach implemented does not work (yet) because of the limitation mentioned in the previous paragraph:

Physics:-`^`(Physics:-Commutator(a, b)+c, 2)

Physics:-`^`(Physics:-Commutator(a, b)+c, 2)

(15)

expand(Physics:-`^`(Physics:-Commutator(a, b)+c, 2))

Physics:-`*`(a, b, a, b)-Physics:-`*`(a, Physics:-`^`(b, 2), a)+Physics:-`*`(a, b, c)-Physics:-`*`(b, Physics:-`^`(a, 2), b)+Physics:-`*`(b, a, b, a)-Physics:-`*`(b, a, c)+Physics:-`*`(c, a, b)-Physics:-`*`(c, b, a)+Physics:-`^`(c, 2)

(16)

In this expression, the 1st, 2nd, 4th and 5th terms have the same operands a, b, a, b and then there are four terms containing the operands a, b, c. We do have an idea of how this could be done too ... :) To be there in one of the next Physics updates.

NULL

NULL


Download Physics[Factor].mw

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

Maple's dsolve numeric can solve delay ODEs and DAEs as of Maple 18. However, if I am not wrong, it cannot solve delay equations with a time dependent history. In this post I show two examples.

Example 1:

y1(t) and y2(t) with time dependent history. Use of piecewise helps this problem to be solved efficiently. Hopefully Maple will add history soon in its capability.

Example 2: 

This is a very a complicated stiff problem from immunology. As of now, I believe only Maple can solve this (other than RADAR5 from Prof. Hairer). Details and plots are posted in the attached code.

 

Let me know if any one has a delay problem that needs to be solved. I have tested many delay problems in Maple (they work fine). The attached examples required addtional tweaking, hence the post.

 

I want to take this opportunity to congratulate and thank Maple's dsolve numeric/delay solvers for their fantastic job. Maple is world leader not because of example1, but because of its ability to solve example 2.

 

 

restart;

 This code is written by Dayaram Sonawane and Venkat R. Subramnian, University of Washington. You will need Maple 18 or later for this. For those who are wanting to solve these problems in earlier versions, I can help them by offering a procedure based approach (less efficient).

Example1 The first example solved is a state dependent delay problem (http://www.mathworks.com/help/matlab/math/state-dependent-delay-problem.html).

 

eq1:= diff(y1(t),t)=y2(t);

eq1 := diff(y1(t), t) = y2(t)

(1)

eq2:=diff(y2(t),t)=-y2(exp(1-y2(t)))*y2(t)^2*exp(1-y2(t));

eq2 := diff(y2(t), t) = -y2(exp(1-y2(t)))*y2(t)^2*exp(1-y2(t))

(2)

 Both y1(t) and y2(t) have time dependent history (y1(t)=log(t) and y2(t)=1/t, t<-0.1). If I am not mistaken one cannot solve this directly using Maple's dsolve numeric command. However, a simple trick can be used to redefine the equations for y1(t) and y2(t) as below

eq3:=diff(y1(t),t)=piecewise(t<=0.1,1/t,y2(t));

eq3 := diff(y1(t), t) = piecewise(t <= .1, 1/t, y2(t))

(3)

eq4:=diff(y2(t),t)=piecewise(t<=0.1,-1/t^2,-y2(exp(1-y2(t)))*y2(t)^2*exp(1-y2(t)));

eq4 := diff(y2(t), t) = piecewise(t <= .1, -1/t^2, -y2(exp(1-y2(t)))*y2(t)^2*exp(1-y2(t)))

(4)

 The problem is solved from a small number close to t = 0 (1e-4) to make Maple's dsolve numeric remember the history till t = 0.1

epsilon:=1e-4;

epsilon := 0.1e-3

(5)

sol:=dsolve({eq3,eq4,y1(epsilon)=log(epsilon),y2(epsilon)=1/epsilon},type=numeric,delaymax=5):

with(plots):

odeplot(sol,[t,y1(t)],0.1..5,thickness=3,axes=boxed);

 

odeplot(sol,[t,y2(t)],0.1..5,thickness=3,axes=boxed);

 

sol(5.0);log(5.0);1/5.0;

[t = 5.0, y1(t) = 1.60942323180838, y2(t) = .199998786891688]

1.609437912

.2000000000

(6)

Tweaking the tolerances and epsilon will get the solution even more closer to the expected answers.

 

 

 Example 2

 The next problem discussed is very stiff, complicated and as of today, according Professor Hairer (one of the world's leading authority in numerical solutions of ODEs, DAEs), cannot be solved by any other code other than his RADAR (5th order implicit Runge Kutta modified for delay equations, Guglielmi N. and Hairer E. (2001) Implementing Radau IIa methods for stiff delay differential equations. Computing 67:1-12). This problem requires very stringent tolerances. For more information read, http://www.scholarpedia.org/article/Stiff_delay_equations. I can safely say that Maple can boast that it can solve this delay differential equation by using a switch function (instead of Heaviside/picecewise function). Code is attached below and results are compared with the output from RADAR code.  Note that dsolve/numeric is probably taking more time steps compared to RADAR, but the fact that Maple's dsolve numeric solved this model (which cannot be solved in Mathematica or MATLAB[needs confirmation for MATLAB]) should make Maple's code writers proud. It is very likely that we will be trying to submit an educational/research article on this topic/example soon to a journal. For some weird reasons, stiff=true gives slightly inaccurate results.

restart:

 

radar5data:=readdata("C:\\Users\\Venkat16core-office\\Google Drive\\waltmanproblem\\sol.txt",[string,string,float,string,string,float,float,float,float,float,float]):

nops(radar5data);

1059

(7)

radar5data[1059];

["X", "=", 300.000000, "Y", "=", 0.6154486288e-15, 0.3377120916e-6, 0.4221403310e-6, 0.2142554563e-5, 299.9999999, 299.6430338]

(8)

eq[1]:=diff(y[1](t),t)=-r*y[1](t)*y[2](t)-s*y[1](t)*y[4](t);

eq[1] := diff(y[1](t), t) = -r*y[1](t)*y[2](t)-s*y[1](t)*y[4](t)

(9)

eq[2]:=diff(y[2](t),t)=-r*y[1](t)*y[2](t)+alpha*r*y[1](y[5](t))*y[2](y[5](t))*H1;#Heaviside(t-35);

eq[2] := diff(y[2](t), t) = -r*y[1](t)*y[2](t)+alpha*r*y[1](y[5](t))*y[2](y[5](t))*H1

(10)

eq[3]:=diff(y[3](t),t)=r*y[1](t)*y[2](t);

eq[3] := diff(y[3](t), t) = r*y[1](t)*y[2](t)

(11)

eq[4]:=diff(y[4](t),t)=-s*y[1](t)*y[4](t)-gamma1*y[4](t)+beta*r*y[1](y[6](t))*y[2](y[6](t))*H2;#Heaviside(t-197);

eq[4] := diff(y[4](t), t) = -s*y[1](t)*y[4](t)-gamma1*y[4](t)+beta*r*y[1](y[6](t))*y[2](y[6](t))*H2

(12)

eq[5]:=diff(y[5](t),t)=H1*(y[1](t)*y[2](t)+y[3](t))/(y[1](y[5](t))*y[2](y[5](t))+y[3](y[5](t)));#eq[7]:=y[7](t)=HH(t);

eq[5] := diff(y[5](t), t) = H1*(y[1](t)*y[2](t)+y[3](t))/(y[1](y[5](t))*y[2](y[5](t))+y[3](y[5](t)))

(13)

eq[6]:=diff(y[6](t),t)=H2*(10.^(-12)*0+y[2](t)+y[3](t))/(10.^(-12)*0+y[2](y[6](t))+y[3](y[6](t)));

eq[6] := diff(y[6](t), t) = H2*(y[2](t)+y[3](t))/(y[2](y[6](t))+y[3](y[6](t)))

(14)

H1:=1/2+1/2*tanh(100*(t-35));H2:=1/2+1/2*tanh(100*(t-197));

H1 := 1/2+(1/2)*tanh(100*t-3500)

H2 := 1/2+(1/2)*tanh(100*t-19700)

(15)

alpha:=1.8;beta:=20.;gamma1:=0.002;r:=5.*10^4;s:=10.^5;

alpha := 1.8

beta := 20.

gamma1 := 0.2e-2

r := 50000.

s := 100000.

(16)

seq(eq[i],i=1..6);

diff(y[1](t), t) = -50000.*y[1](t)*y[2](t)-100000.*y[1](t)*y[4](t), diff(y[2](t), t) = -50000.*y[1](t)*y[2](t)+90000.0*y[1](y[5](t))*y[2](y[5](t))*(1/2+(1/2)*tanh(100*t-3500)), diff(y[3](t), t) = 50000.*y[1](t)*y[2](t), diff(y[4](t), t) = -100000.*y[1](t)*y[4](t)-0.2e-2*y[4](t)+1000000.*y[1](y[6](t))*y[2](y[6](t))*(1/2+(1/2)*tanh(100*t-19700)), diff(y[5](t), t) = (1/2+(1/2)*tanh(100*t-3500))*(y[1](t)*y[2](t)+y[3](t))/(y[1](y[5](t))*y[2](y[5](t))+y[3](y[5](t))), diff(y[6](t), t) = (1/2+(1/2)*tanh(100*t-19700))*(y[2](t)+y[3](t))/(y[2](y[6](t))+y[3](y[6](t)))

(17)

ics:=y[1](0)=5.*10^(-6),y[2](0)=10.^(-15),y[3](0)=0,y[4](0)=0,y[5](0)=1e-40,y[6](0)=1e-20;

ics := y[1](0) = 0.5000000000e-5, y[2](0) = 0.1000000000e-14, y[3](0) = 0, y[4](0) = 0, y[5](0) = 0.1e-39, y[6](0) = 0.1e-19

(18)

#infolevel[all]:=10;

sol:=dsolve({seq(eq[i],i=1..6),ics},type=numeric,delaymax=300,initstep=1e-6,abserr=[1e-21,1e-21,1e-21,1e-21,1e-9,1e-9],[y[1](t),y[2](t),y[3](t),y[4](t),y[5](t),y[6](t)],relerr=1e-9,maxstep=10,optimize=false,compile=true,maxfun=0):

 

 

 note that compile = true was used for efficiency

t11:=time():sol(300);time()-t11;

[t = 300., y[1](t) = 0.615611371327094e-15, y[2](t) = 0.337706811581908e-6, y[3](t) = 0.422136411682798e-6, y[4](t) = 0.214253771204037e-5, y[5](t) = 299.999986716780, y[6](t) = 299.643054284209]

.141

(19)

with(plots):

nd:=nops(radar5data);

nd := 1059

(20)

radar5data[nd];

["X", "=", 300.000000, "Y", "=", 0.6154486288e-15, 0.3377120916e-6, 0.4221403310e-6, 0.2142554563e-5, 299.9999999, 299.6430338]

(21)

 Values at t = 300 match with expected results.

pr[1]:=plot([seq([radar5data[i][3],log(radar5data[i][6])/log(10)],i=1..nd)],style=point,color=green):

p[1]:=odeplot(sol,[t,log(y[1](t))/log(10)],0..300,axes=boxed,thickness=3):

display({pr[1],p[1]});

 

pr[2]:=plot([seq([radar5data[i][3],log(radar5data[i][7])/log(10)],i=1..nd)],style=point,color=green):

p[2]:=odeplot(sol,[t,log(y[2](t))/log(10)],0..300,axes=boxed,thickness=3,numpoints=1000):

display({pr[2],p[2]});

 

pr[3]:=plot([seq([radar5data[i][3],log(radar5data[i][8])/log(10)],i=2..nd)],style=point,color=green):

 

p[3]:=odeplot(sol,[t,log(y[3](t))/log(10)],0..300,axes=boxed,thickness=3):

display({pr[3],p[3]});

 

pr[4]:=plot([seq([radar5data[i][3],log(radar5data[i][9])/log(10)],i=496..nd)],style=point,color=green,view=[197..300,-9..-5]):

p[4]:=odeplot(sol,[t,log(y[4](t))/log(10)],197..300,axes=boxed,thickness=3,view=[197..300,-9..-5]):

display({pr[4],p[4]});

 

pr[5]:=plot([seq([radar5data[i][3],radar5data[i][10]],i=1..nd)],style=point,color=green):

p[5]:=odeplot(sol,[t,y[5](t)],0..300,axes=boxed,thickness=3):

display({pr[5],p[5]});

 

pr[6]:=plot([seq([radar5data[i][3],radar5data[i][11]],i=1..nd)],style=point,color=green):

p[6]:=odeplot(sol,[t,y[6](t)],0..300,axes=boxed,thickness=3):

display({pr[6],p[6]});

 


Download delayimmunetopost.mws

As a very good application for viewing and calculation of the components of acceleration either tangential or normal. Besides immediately it is shown an Application for physics.

Componentes_de_la_Acelelación.mw

(in spanish)

L. Araujo C.

 

 

Here we see the projection of a vector onto another using different concepts ranging from linear algebra to vector calculus. Implemented components thus seen in three-dimensional space.

 

Proyecciones_Vectoriales.mw

(in spanish)

L.Araujo C.

Here the potential of maple 2015 to the quantitative study of the decomposition of a vector table is shown in two dimensions. Application for the exclusive use of engineering students, which was implemented with embedded components.

Atte.

Lenin Araujo Castillo

Archivo Corregido:  Decomposición_Vectorial_Corregido.mw

The first instalment in the Hollywood Math webinar series will be returning live this September 17th! It will be presented by Daniel Skoog, our Maple Product Manager.

Over its storied and intriguing history, Hollywood has entertained us with many mathematical moments in film. John Nash in “A Beautiful Mind,” the brilliant janitor in “Good Will Hunting,” the number theory genius in “Pi,” and even Abbott and Costello are just a few of the Hollywood “mathematicians” that come to mind.

During this webinar we will present a number of examples of mathematics in film, including those done capably, as well as questionable and downright “creative” treatments. See relevant, exciting examples that you can use to engage your students, or attend this webinar simply for its entertainment value. Have you ever wondered if the bus could really have jumped the gap in “Speed?” We’ve got the answer! Anyone with an interest in mathematics, especially high school and early college math educators, will be both entertained and informed by attending this webinar.

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


In connection with recent developments for symbolic sequences, a number of improvements were implemented regarding symbolic differentiation, that is the computation of n^th order derivatives were n is a symbol, the simplest example being the n^th derivative of the exponential, which of course is the exponential itself. This post is about these developments, done in collaboration with Katherina von Bülow, and available for download as usual from the Maplesoft R&D web page for Differential Equations and Mathematical functions (the update itself is bundled with the official updates of the Maple Physics package).

 

It is important to note that Maple is pioneer in having an actual implementation of symbolic differentiation, something that works for real, since several releases.  The development, however, was somewhat stuck because we were unable to compute the symbolic n^th derivative of a composite function f(g(z)). A formula for this problem is actually known, it is the Faà di Bruno formula, but, in order to implement it, first we were missing the incomplete Bell functions , that got implemented in Maple 15, nice, but then we were still missing differentiating symbolic sequences, and functions whose arguments are symbolic sequences (i.e. the number of arguments of the function is n, a symbol, of unknown value at the time of differentiating). All this got implemented now within the new MathematicalFunctions:-Sequence package, opening the door widely to these improvements in n^th differentiation.

 

The symbolic differentiation code works as mostly all other computer algebra code, by mapping complicated problems into a composition of simpler problems all of which are tractable; what follows is then an illustration of these basic cases.

 

Among the simplest new case that can now be handled there is that of a power where the exponent is linear in the differentiation variable. This is actually an easy problem

(%diff = diff)(f^(alpha*z+beta), `$`(z, n))

%diff(f^(alpha*z+beta), `$`(z, n)) = alpha^n*f^(alpha*z+beta)*ln(f)^n

(1)

More complicated, consider the k^th power of a generic function; the corresponding symbolic derivative can be mapped into a sum of symbolic derivatives of powers of g(z) with lower degree

(%diff = diff)(g(z)^k, `$`(z, n))

%diff(g(z)^k, `$`(z, n)) = k*binomial(n-k, n)*(Sum((-1)^_k1*binomial(n, _k1)*g(z)^(k-_k1)*(Diff(g(z)^_k1, [`$`(z, n)]))/(k-_k1), _k1 = 0 .. n))

(2)

In some cases where g(z) is a known function, the computation can be carried on furthermore. For example, for g = ln the result can be expressed using Stirling numbers of the first kind

(%diff = diff)(ln(alpha*z+beta)^k, `$`(z, n))

%diff(ln(alpha*z+beta)^k, `$`(z, n)) = alpha^n*(Sum(pochhammer(k-_k1+1, _k1)*Stirling1(n, _k1)*ln(alpha*z+beta)^(k-_k1), _k1 = 0 .. n))/(alpha*z+beta)^n

(3)

The case of sin and cos are relatively simpler, but then assumptions on the exponent are required in order to proceed further ahead from (2), for example

`assuming`([(%diff = diff)(sin(alpha*z+beta)^k, `$`(z, n))], [k::posint])

%diff(sin(alpha*z+beta)^k, `$`(z, n)) = (-1)^k*piecewise(n = 0, (-sin(alpha*z+beta))^k, alpha^n*I^n*(Sum(binomial(k, _k1)*(2*_k1-k)^n*exp(I*(2*_k1-k)*(alpha*z+beta+(1/2)*Pi)), _k1 = 0 .. k))/2^k)

(4)

The case of functions of arbitrary number of variables (typical situation where symbolic sequences are required) is now handled properly. This is the pFq hypergeometric function of symbolic order p and q 

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

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

(5)

The case of the MeijerG function is more complicated, but in practice, for the computer, once it knows how to handle symbolic sequences, the more involved problem becomes computable

(%diff = diff)(MeijerG([[`$`(a[i], i = 1 .. n)], [`$`(b[i], i = n+1 .. p)]], [[`$`(b[i], i = 1 .. m)], [`$`(b[i], i = m+1 .. q)]], z), `$`(z, k))

%diff(MeijerG([[`$`(a[i], i = 1 .. n)], [`$`(b[i], i = n+1 .. p)]], [[`$`(b[i], i = 1 .. m)], [`$`(b[i], i = m+1 .. q)]], z), `$`(z, k)) = MeijerG([[-k, `$`(a[i]-k, i = 1 .. n)], [`$`(b[i]-k, i = n+1 .. p)]], [[`$`(b[i]-k, i = 1 .. m)], [0, `$`(b[i]-k, i = m+1 .. q)]], z)

(6)

Not only the mathematics of this result is correct: the object returned is actually computable to the end (if you provide the values of n, p, m and q), and the typesetting is actually fully readable, as in textbooks, including copy and paste working properly; all this is new.

The n^th derivative of a number of mathematical functions that were not implemented before, are now also implemented, covering the gaps, for example:

(%diff = diff)(BellB(a, z), `$`(z, n))

%diff(BellB(a, z), `$`(z, n)) = Sum(Stirling2(a, _k1)*pochhammer(_k1-n+1, n)*z^(_k1-n), _k1 = 0 .. a)

(7)

(%diff = diff)(bernoulli(z), `$`(z, n))

%diff(bernoulli(z), `$`(z, n)) = pochhammer(nu-n+1, n)*bernoulli(nu-n, z)

(8)

(%diff = diff)(binomial(z, m), `$`(z, n))

%diff(binomial(z, m), `$`(z, n)) = (Sum((-1)^(_k1+m)*Stirling1(m, _k1)*pochhammer(_k1-n+1, n)*(z-m+1)^(_k1-n), _k1 = 1 .. m))/factorial(m)

(9)

(%diff = diff)(euler(a, z), `$`(z, n))

%diff(euler(a, z), `$`(z, n)) = pochhammer(a-n+1, n)*euler(a-n, z)

(10)

In the same way the fundamental formulas for the n^th derivative of all the 12 elliptic Jacobi functions  as well as the four elliptic JacobiTheta functions,  the LambertW , LegendreP  and some others are now all implemented.

Finally there is the "holy grail" of this problem: the n^th derivative of a composite function f(g(z)) - this always-unreachable implementation of Faa di Bruno formula. We now have it :)

(%diff = diff)(f(g(z)), `$`(z, n))

%diff(f(g(z)), `$`(z, n)) = Sum(((D@@k)(f))(g(z))*IncompleteBellB(n, k, `$`(diff(g(z), [`$`(z, j)]), j = 1 .. n-k+1)), k = 0 .. n)

(11)

Note the symbolic sequence of symbolic order derivatives of lower degree, both of of f and g, also within the arguments of the IncompleteBellB function. This is a very abstract formula ... And does this really work? Of course it does :). Consider, for instance, a case where the n^th derivatives of f(z) and g(z) can both be computed by the system:

sin(cos(alpha*z+beta))

sin(cos(alpha*z+beta))

(12)

This is the n^th derivative expressed using Faa di Bruno's formula, in turn expressed using symbolic sequences within the IncompleteBellB  function

(%diff = diff)(sin(cos(alpha*z+beta)), `$`(z, n))

%diff(sin(cos(alpha*z+beta)), `$`(z, n)) = Sum(sin(cos(alpha*z+beta)+(1/2)*k*Pi)*IncompleteBellB(n, k, `$`(cos(alpha*z+beta+(1/2)*j*Pi)*alpha^j, j = 1 .. n-k+1)), k = 0 .. n)

(13)

These results can all be verified. Take for instance n = 3

eval(%diff(sin(cos(alpha*z+beta)), `$`(z, n)) = Sum(sin(cos(alpha*z+beta)+(1/2)*k*Pi)*IncompleteBellB(n, k, `$`(cos(alpha*z+beta+(1/2)*j*Pi)*alpha^j, j = 1 .. n-k+1)), k = 0 .. n), n = 3)

%diff(sin(cos(alpha*z+beta)), z, z, z) = Sum(sin(cos(alpha*z+beta)+(1/2)*k*Pi)*IncompleteBellB(3, k, `$`(cos(alpha*z+beta+(1/2)*j*Pi)*alpha^j, j = 1 .. 4-k)), k = 0 .. 3)

(14)

Compute now the inert functions: on the left-hand side this is just the (now explicit) 3rd order derivative, while on the right-hand side we have a sum of IncompleteBellB  functions, where the number of arguments, expressed in (13) using symbolic sequences that depend on the summation index k and the differentiation order n, now in (14) depend only on k, and get transformed into explicit sequences of arguments when the summation is performed and k assumes integer values

value(%diff(sin(cos(alpha*z+beta)), z, z, z) = Sum(sin(cos(alpha*z+beta)+(1/2)*k*Pi)*IncompleteBellB(3, k, `$`(cos(alpha*z+beta+(1/2)*j*Pi)*alpha^j, j = 1 .. 4-k)), k = 0 .. 3))

alpha^3*sin(alpha*z+beta)*cos(cos(alpha*z+beta))-3*alpha^3*cos(alpha*z+beta)*sin(alpha*z+beta)*sin(cos(alpha*z+beta))+alpha^3*sin(alpha*z+beta)^3*cos(cos(alpha*z+beta)) = alpha^3*sin(alpha*z+beta)*cos(cos(alpha*z+beta))-3*alpha^3*cos(alpha*z+beta)*sin(alpha*z+beta)*sin(cos(alpha*z+beta))+alpha^3*sin(alpha*z+beta)^3*cos(cos(alpha*z+beta))

(15)

Take left-hand side minus right-hand side

simplify((lhs-rhs)(alpha^3*sin(alpha*z+beta)*cos(cos(alpha*z+beta))-3*alpha^3*cos(alpha*z+beta)*sin(alpha*z+beta)*sin(cos(alpha*z+beta))+alpha^3*sin(alpha*z+beta)^3*cos(cos(alpha*z+beta)) = alpha^3*sin(alpha*z+beta)*cos(cos(alpha*z+beta))-3*alpha^3*cos(alpha*z+beta)*sin(alpha*z+beta)*sin(cos(alpha*z+beta))+alpha^3*sin(alpha*z+beta)^3*cos(cos(alpha*z+beta))))

0

(16)

NULL

:)


Download SymbolicOrderDifferentiation.mw


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

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.

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

First 30 31 32 33 34 35 36 Last Page 32 of 59