Maple-MATLAB® Link Application
Initiate the MATLAB® link. Use the with command to access the functions in some useful packages by their short names:
restart: with(LinearAlgebra): with(plots): with(Matlab):
For more information on the Maple-MATLAB® link, see Matlab.
Structural Analysis: A First Approximation
This formulation is used to compute the lowest natural frequencies and modes of a highly idealized 22-story building.
The model equations may be formulated with the following matrix assignments.
The mass matrix M is a matrix with m on the diagonal:
m := 5000.: M := DiagonalMatrix([seq(m,i=1..22)], outputoptions=[datatype=float,shape=symmetric,attributes=[positive_definite],storage=rectangular]);
To examine any of these Matrices, the Structured Data Browser can be used, by right-clicking the output Matrix and selecting Browse.
The Stiffness matrix K is tridiagonal with 2k on the center diagonal, and -k on the adjacent diagonals:
k := 1.25e8: K := BandMatrix([-k,2*k,-k],1,22,22, outputoptions=[datatype=float,shape=symmetric,storage=rectangular]): K[22,22] := -k: K;
The Eigenvalues and Eigenvectors computed with MATLAB® are found:
(P,W):=eig(K,M,eigenvectors=true):
The Eigenvalues are:
W;
The Eigenvectors are:
P;
The equivalent computations in the Maple environment:
(W2,P2) := Eigenvectors(K,M):
W2;
the Eigenvectors are (in no particular order):
P2;
Heat Transfer: Finite Difference Solution
A difference equation method is used to find the static temperature distribution in a flat rectangular plate, given its boundary is held at a fixed temperature profile.
We model the plate as a 3 x 7 grid of nodes, where the nodes may be thought of as being interconnected with a square mesh of heat conductors.
With our 21 plate-internal nodes and 20 boundary conditions ( = 7+3+7+3 ), the finite-difference 2-dimensional Laplace equation gives us an inhomogeneous linear system in 21 unknowns (the internal nodal temperatures). Representing the internal nodal temperatures U[i,j] by the column vector [ U[1,1], ..., U[1,7], U[2,1], ..., U[2,7], U[3,1], ..., U[3,7] ],we have the system AU=B, where the matrix A and column vector B encode the interconnectivity of the nodes and the profile of the boundary temperature:
A := BandMatrix([1,0,0,0,0,0,1,-4,1,0,0,0,0,0,1],7,21,21,outputoptions=[datatype=float[8]]): A[7,8]:=0:A[8,7]:=0: A[14,15]:=0:A[15,14]:=0: A;
T := 100.: B := Vector[column]( 21, [-T,-T,-T,-T,-T,-T,-T,0,0,0,0,0,0,0,0,0,0,0,0,0,0], datatype=float[8] );
We have fixed the temperature along the U[1,1], ..., U[1,7] edge at 100 units, and the other three edges at 0 temperature units.
Now, we solve the system in the MATLAB® environment:
setvar("B", B);
setvar("A", A);
evalM("U=inv(A) * B"); U1 := getvar("U");
Or, we can solve the system entirely within Maple:
U2 := MatrixInverse(A) . B;
Now, let us take the Maple solution vector and re-arrange its values in a 3 x 7 matrix that corresponds to the actual rectangular plate:
P := <<U2[1..7]>|<U2[8..14]>|<U2[15..21]>>:
plots[matrixplot](P,heights=histogram,axes=boxed,labels=["","","T"],title="temperature distribution");
Simulation of a Linear Oscillator
A simple spring-mass-dashpot is modeled as a second-order linear oscillator. The simulation equations are coded as MATLAB® function files and are called from the Maple environment by using the ode45 command.
The MATLAB® function stored in the file mass_eqn.m can be created within Maple as follows, or it can be created by using any text editor. Here, we create the file in the current directory.
currentdir(kernelopts(homedir));
file := open("mass_eqn.m", WRITE):
writeline(file, "function xdot=mass_eqn(t,x)"): writeline(file, " global M C K"): writeline(file, " xdot = [ x(2); (1/M*(-K*x(1)-C*x(2))) ];"):
writeline(file, "end"):
close(file);
The oscillator parameters of Mass M, Damping C, and Stiffness K are:
M:=100; C:=10; K:=8;
M:=100
C:=10
K:=8
Pass the variables to MATLAB®. We must enforce the above variables as global in the MATLAB® environment so that the function mass_eqn can use them.
setvar("M", M, 'globalvar'): setvar("C", C, 'globalvar'): setvar("K", K, 'globalvar'):
Allow MATLAB® to solve the ODE.
(t,x) := ode45("mass_eqn", 0..50, [2,0]);
plots[pointplot]([seq([t[i],x[i,1]],i=1..Dimensions(t))], symbol=DIAMOND, symbolsize=8);
Data Analysis: Fast Fourier Transform (FFT)
A set of experimental numerical data can be analyzed with MATLAB® as a numerical engine.
t := [seq(i*.001, i=0..600)]:
tmax := nops(t);
tmax:=601
x := map(proc(T) evalhf(sin(2*Pi*50*T) + sin(2*Pi*120*T)) end, t):
data := map(proc(X) X+2*stats[random,normald](1) end, x):
Y := fft(data,512);
Pyy := [seq(abs(Y[i])^2/ 512, i=1..512)]:
freq := [seq(1000*i/512., i=0..511)]:
plots[pointplot]([seq([t[i],data[i]],i=1..nops(t))],symbol=DIAMOND,symbolsize=7);
plots[pointplot]([seq([freq[i],Pyy[i]],i=1..512)],style=line);
Alternatively, the same calculations can be done by using MATLAB® syntax almost entirely.
evalM("t=0:0.001:0.6");
t:=getvar("t");
tmax := op(2,size(t));
evalM("x=sin(2*pi*50*t) + sin(2*pi*120*t)");
evalM("y=x+2*randn(size(t))");
data:=getvar("y");
evalM("Y=fft(y,512)");
evalM("Pyy=Y.*conj(Y)/512");
Pyy:=getvar("Pyy");
evalM("f=1000*(0:511)/512");
freq:=getvar("f"):
Dimensions(data);
601
Dimensions(freq);
512
Dimensions(Pyy);
plots[pointplot]([seq([t[i],data[i]],i=1..tmax)],symbol=DIAMOND,symbolsize=7);
Return to Index for Example Worksheets
Download Help Document
