Scot Gould

Scot Gould

872 Reputation

14 Badges

11 years, 223 days
Claremont McKenna, Pitzer, Scripps College
Professor of Physics
Upland, California, United States
Dr. Scot Gould is a professor of physics in the W.M. Keck Science Department of Claremont McKenna, Pitzer, and Scripps Colleges - members of The Claremont Colleges in California. He was involved in the early development of the atomic force microscope. His research has included numerous studies and experiments using scanning probe microscopes, particularly those involving natural fibers such as spider silk. More recently, he was involved in developing and sustaining AISS. This full-year multi-unit, non-traditional, interdisciplinary undergraduate science education course integrated topics from biology, chemistry, physics, mathematics, and computer science. His current interest is integrating computational topics into the physics curriculum. He teaches the use of Maple's computer algebraic and numerical systems to assist students in modeling and visualizing physical and biological systems. His Dirac-notation-based quantum mechanics course is taught solely through Maple.

MaplePrimes Activity


These are questions asked by Scot Gould

dsolve gets stuck on a problem. Attached is the worksheet.

In Case I, dsolve solves the two ODEs separately. In Case II, dsolve solves them together, but only if I define a constant. In Case III the constant is not defined. dsolve never returns, and I have to hit the stop sign. The question I have: Is there a consistent method to the failure so that I can avoid the problem in the future with other ODEs?

dsolve_bug.mw

Question: Is there a style for the icons from the palettes? 

(The rest of this post is a screenshot because MaplePrimes can't show the change in fonts.)s a

The following is not a profound problem, and there is an obvious solution,

but it came up, and I would like to learn more about it.

 

Even though I recommend the add procedure when summing up individual entities,

my students keep showing me how smart the sum procedure is. Which makes

our worksheets more readable and reproducible for Maple users who are less frequent.

 

For example:

 

restart; Xlist := [1, 2, 3]; N := numelems(Xlist)

3

Using palette icon:

sum(Xlist[n], n = 1 .. N)

6

Cool!  Which means

sum(Xlist[n], n = 1 .. N)

6

But if we use the same palette icon for a vector

Xvector := convert(Xlist, Vector); sum(Xvector[n], n = 1 .. N)

Error, bad index into Vector

Because I believe this fails

sum(Xvector[n], n = 1 .. N)

Error, bad index into Vector

 

Would someone please teach me how I can see why the sum of a list

works, but does the sum of a vector fail?

Download MaplePrimes_sum_list_vector.mw

Since my problem involving a 2d heat equation was answered so well,

(hats off to nm !!) I thought I would try another PDE question.

 

Question: can pdsolve handle a combination of the function and its derivative for

the boundary conditions, i.e., Robin boundary conditions?

 

I looked through the document example, pdsolve_boundaryconditions, which was

beautifully written. However, I saw no example with Robin BCs.

 

"restart; interface(imaginaryunit = i):      StartTemp(x):= `T__0`*(e)^(-`alpha__0` )*cos(Pi*(x)/(`L__0`));   `L__0` := 10: #` Length of rod` `T__0` := 100: #` max temperature` `alpha__0` := 6/(100): #`  Dampening factor` plot(StartTemp(x) , x = 0 .. `L__0`,  thickness = 5,             view = [0..`L__0`, -100..100],  labels = ["Position on rod", "Temperature"], labelfont = [Times, 14],             labeldirections = [horizontal, vertical], title = "Temperature in rod at t = 0" , titlefont = [Times, 16],             size = [500, 225]);"

proc (x) options operator, arrow, function_assign; T__0*exp(-alpha__0)*cos(Pi*x/L__0) end proc

 

 

The equation and the initial condition:

heat_equation := diff(T(x, t), t) = k*(diff(T(x, t), x, x)); initial_condition := T(x, 0) = T__i(x)

Using a linear combination of the function and its derivative, i.e., Robin BCs.

boundary_conditions := `α__r`*T(0, t)+`β__r`*(D[1](T))(0, t) = f(t), `α__r`*T(L__0, t)+`β__r`*(D[1](T))(L__0, t) = g(t)

 

Setting the constants

"  `alpha__r` := 1:  #` Constant for boundary condition`  `beta__r` := 2: #` Constant for boundary condition`     `T__s`:= 10:  #` Time interval`  k:=9/(10):  #` Heat-conductivity of material`    f(t) :=-1/(10) t:                  g(t):=1/(10) t:  "

Solving using nm's technique of not defining the function before calling pdsolve

sols := pdsolve({heat_equation, boundary_conditions, initial_condition}, T(x, t))

 

Now include the initial temperature function

"`T__i`(x) := StartTemp(x):  "

Extracting the solution, and using 20 terms

sol := eval(rhs(sols), infinity = 20)

-134/135+(1/270)*x^3-(7/90)*x^2+(67/135)*x+Sum(exp(-(9/1000)*Pi^2*n^2*t)*(Pi*n*cos((1/10)*n*Pi*x)-5*sin((1/10)*n*Pi*x))*(Int(-(x^3-21*x^2+134*x-27000*exp(-3/50)*cos((1/10)*Pi*x)-268)*(Pi*n*cos((1/10)*n*Pi*x)-5*sin((1/10)*n*Pi*x)), x = 0 .. 10))/(1350*Pi^2*n^2+33750), n = 1 .. 20)+(1/50)*t*x-(7/50)*t

(1)

Evaluate the integrals and build the solution to the temperature as a function of position and time

T__sol := MakeFunction(value(sol), x, t)

 

Here is the problem:

in plotting the solution at t = 0, we do NOT reproduce the initial temperature function.

 

plot(T__sol(x, 0), x = 0 .. L__0, -120 .. 120, thickness = 5, size = [500, 225])

 

If I use boundary conditions for either Dirichlet or Neumann, i.e., T(0, t) & T(L__0, t) 

or "`T__x`(0,t) & `T__x`(`L__0`,t)", the initial function is reproduced at t = 0.

 

Is this example one of those situations that was addressed in the notes?

Given a problem with PDEs and some BCs, note the following:

• 

Depending on how you write a general PDE solution, it becomes
possible or nearly impossible to adjust it to match the BCs.

 

NULL

Download Heat_Equation_1D_Robin_Conditions.mw

2D Heat Equation problem.

 

In most cases, pdsolve generates expected results. For this one, it does not. Is there an obvious error

in my worksheet?

 

restart; interface(imaginaryunit = i); heat_equation := diff(T(x, y, t), t) = k*(diff(T(x, y, t), x, x)+diff(T(x, y, t), y, y))

 

Writing out the boundary conditions:

bces := T(0, y, t) = 0, T(L, y, t) = 0, T(x, 0, t) = 0, T(x, L, t) = 0

Writing out the initial conditions

ice := T(x, y, 0) = T__0(x, y)

 

where

"`T__0`(x,y) := 1/(4)*sin((2 Pi*x)/(L) )^(2)*sin((3 Pi*y)/(L)):"

 

Showing the initial condition of the system.

L := 1; plot3d_display := scaling = constrained, size = [600, 600], orientation = [-121, 75, 1]; plot3d(T__0(x, y), x = 0 .. L, y = 0 .. L, plot3d_display)

 

Solving the heat equation.

sols := `assuming`([pdsolve({bces, ice, heat_equation}, T(x, y, t))], [k > 0])

T(x, y, t) = Sum(Sum(piecewise(n = 4, -((1/8)*I)*sin(n1*Pi*y)*sin(4*Pi*x)*exp(-Pi^2*k*t*(n1^2+16)), -4*sin(n*Pi*x)*exp(-Pi^2*k*t*(n^2+n1^2))*(-1+(-1)^n)*sin(n1*Pi*y)/(Pi*(n^3-16*n))), n = 1 .. infinity), n1 = 1 .. infinity)

(1)

Notice the imaginary expression in the solution?

 

Tsol := value(eval(rhs(sols), infinity = 12))

 

Looking at one value that includes an imaginary term.

eval(Tsol, {k = .1, t = 0, x = .2, y = .2})

-.3490991334-.1130635622*I

(2)

And plotting at t = 0.

plot3d(eval(Tsol, {k = .3, t = 0}), x = 0 .. L, y = 0 .. L, plot3d_display)

 

NULL

Note, if the sin(3*Pi*y/L) term is squared, then pdsolve returns a real solution.

NULL

Download Heat_Eq_2D.mw

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