Carl Love

Carl Love

28035 Reputation

25 Badges

12 years, 319 days
Himself
Wayland, Massachusetts, United States
My name was formerly Carl Devore.

MaplePrimes Activity


These are answers submitted by Carl Love

The corresponding input is Eval(H1, t= 0.45). That's with a capital E so that it's inert. To activate the evaluation, you'd use value(%) (if you were actually running Maple, of course).

Use the -q (for quiet) flag on the cmaple command, for example

cmaple -q < myinputfile.mpl

In your module, include a procedure named ModuleLoad (you must use precisely this name). It may be declared local or export, but it needs to be one of these. Inside this procedure, include the statement :-Clr:= B. Right before the end module, include the statement ModuleLoad().

I don't understand why you want Clr to be global rather than a module export.

You need to use eval rather than subs. It is a lot like subs, but much more mathematically sophisticated, and it understands how to handle derivatives, integrals, etc. In constrat, subs is intended primarily for low-level programming. It just blindly substitutes what you tell it to.

Here's how to use eval in your situation:

g:= (x,y)-> diff(u1(x, y), x, x) + diff(u2(x, y), x, y):
f:= (y)-> eval(g(x,Y), [x= 0, Y= y]);

Note the distinction that I have made between lowercase y and uppercase Y.

Here's a completely different way:

g:= D[1,1](u1) + D[1,2](u2):
f:= (y)-> g(0,y):

 

 

 

Here's an elegant and automated solution to your problem. The main idea is to create a new ODE and a new BC by differentiating the expression that you want evaluated.
 

An iterative solution to a parametrized BVP

Carl J Love 2016-Dec-22


restart:


#Your original system:
#(This should be the only thing you that need to modify):

EQs:=
   diff(f(x),x$3) + f(x)*diff(f(x),x$2) - a*diff(f(x),x$1)^2 = 0,
   diff(g(x),x$2) + b*f(x)*diff(g(x),x$1) = 0
;
LB:= 0:  UB:= 5: #upper and lower boundaries
BCs:= f(LB)=0, D(f)(LB)=1, D(f)(UB)=0, g(LB)=0.5, g(UB)=0;
expr:= a*f(x) + b*diff(g(x),x) + c*diff(f(x),x)*g(x);
#expr to be evaluated at:
x__e:= 1.; Params:= [a,b,c]=~ ([seq(0..1, 0.2)] $ 3);

diff(diff(diff(f(x), x), x), x)+f(x)*(diff(diff(f(x), x), x))-a*(diff(f(x), x))^2 = 0, diff(diff(g(x), x), x)+b*f(x)*(diff(g(x), x)) = 0

f(0) = 0, (D(f))(0) = 1, (D(f))(5) = 0, g(0) = .5, g(5) = 0

a*f(x)+b*(diff(g(x), x))+c*(diff(f(x), x))*g(x)

1.

[a = [0, .2, .4, .6, .8, 1.0], b = [0, .2, .4, .6, .8, 1.0], c = [0, .2, .4, .6, .8, 1.0]]


#Implicity define a new function h(x) to be expr, and make it into
#a new ODE:
Expr:= h(x) = expr:
newEQ:= diff(h(x) = Expr, x);

#Create a BC for h(0) by evaluating Expr at x=0 and then evaluating
#wrt BCs:
#(This step can be tricky; may be difficult to fully automate.)
#(Might also want to try to make a new BC at UB if the one at LB
#doesn't work.)
newBC:= eval(convert(eval(Expr, {x= LB}), D), {BCs});

#Abstract the BVP system parametrized by (a,b,c):
sys:= unapply({EQs, newEQ, BCs, newBC}, [a,b,c]):

#And use that to create a parametrized solver that explicitly
#evaluates h(x) at a point.
_H:= (a::realcons, b::realcons, c::realcons, x0::realcons)->
   eval(h('x'), dsolve(sys(a,b,c), numeric)(x0))
:

#Example usage:
h(x__e) = _H(0.5, 0.5, 0.5, x__e);
 

diff(h(x), x) = (diff(h(x), x) = a*(diff(f(x), x))+b*(diff(diff(g(x), x), x))+c*(diff(diff(f(x), x), x))*g(x)+c*(diff(f(x), x))*(diff(g(x), x)))

h(0) = b*(D(g))(0)+.5*c

h(1.) = HFloat(0.31148732145735375)

#Do calculations, putting them into a table:
Calc:= proc()
local abc, ABC, Htable:= table();

   for abc in Iterator:-CartesianProduct(eval([a,b,c], Params)[]) do
      ABC:= convert(abc, list);
      assign(`?[]`(Htable, ABC) = _H(ABC[], x__e))
   end do;
   Htable
end proc:

#Measure resource consumption for the "hard" computation:
Htable:= CodeTools:-Usage(Calc()):

memory used=431.50MiB, alloc change=0 bytes, cpu time=5.44s, real time=5.46s, gc time=218.75ms


#Put the table in an elegantly displayed format:
N:= mul(nops~(eval([a,b,c], Params))): #number of parameter tuples.
sav:= interface(rtablesize= N+2):
DataFrame(
   Matrix((`[]`@op)~([entries(Htable, pairs, indexorder)])),
   columns= [a, b, c, h(x__e)], rows= ['convert(``, `local`)' $ N]
);
interface(rtablesize= sav):

RTABLE(18446744074183065598, anything, Matrix, rectangular, Fortran_order, [], 2, 1 .. 217, 1 .. 5)

 


 

Download ParamtrizedBVP.mw

A direct translation into Maple of what you have written is

local D:= n->
   Matrix(
      ((n+1)$2),
      (i,k)->
         if i <> k then P[n](t[i])/P[n](t[k])/(t[i] - t[k])
         elif i = 1 then -n*(n+1)/2
         elif i = N+1 then n*(n+1)/2
         else 0
         end if
   )
:

I can't read your denominators on n*(n+1). I have guessed 2.

Your method of presentation shows disrespect for this forum and its members. If you can't come up with a better presentation, you should at least first ask for permission to post something so crude, and beg for the forgiveness of your readers.

When you make an assignment to F[1], it also changes the meaning of F, which becomes a table, not a vector field. If you want a new variable which displays as F with a subscript, use F__1 (that's two underscores). That is a completely new variable, with no connection to F. When any name is followed by [anything], the resulting structure has some relation to the parent name.

Regarding your title: Consider this simple example:

x:= 2:
x:= 1/(x-1);
x:= 1/(x-1);

Surely you can understand that it is correct that the second command works and that the third one returns an error although they are syntacticly the same command.

One way to do what you want (and it's certainly the easiest programmatic way) is make a Matrix whose first row is the headers and whose second row is the matrices. So, some of the elements of this Matrix are themselves matrices.

Display:= Matrix(2,3, [[`Matrix 1`, `Matrix 2`, `Matrix 3`]]):
Display[2,1]:= A1: Display[2,2]:= A2: Display[2,3]:= A3:
Display;

    

If A and B are vectors, matrices, arrays, etc., of the same size and shape, then A *~ B performs the elementwise multiplication.

How would you recommend that this feature be indexed in the help pages, since you seem to have tried to find it? It is listed under ?elementwise

In Maple, e^x is exp(x). There is no caret (^) in it.

See ?Student,Calculus1,VolumeOfRevolution and ?Student,Calculus1,VolumeOfRevolutionTutor

A surface integral is a method of computing the flux (or flow) of a vector field (such as a force field or a moving fluid) through a surface. It is not ordinarily used to compute volume.

If you parametrize a surface, then in any integral that uses that parametrization, you need to multiply the integrand by the determinant of the Jacobian (the matrix of partial derivatives) of the transformation. This is completely analogous to changing dx to du when doing a u-substitution in a single integral. This determinant is where the extra factor of r comes from for polar and cylindrical coordinates and the extra factor of rho^2*sin(theta) for spherical coordinates. You haven't done this in your integral. I'm not saying that you should go do it! Because the easiest way to do this volume integral is in cylindrical coordinates.

Your plot is pretty good. It definitely shows the correct two surfaces.

In the worksheet below, I compute the volume five different ways, all using integrals. The Answer by Tom Leslie is incorrect.
 

restart:

plots:-display(
   plot3d( #upper hemisphere
      [r, theta, sqrt(4-r^2)], r= 1..2, theta= 0..2*Pi,
      coords= cylindrical, style= wireframe, color= brown, thickness= 3
   ),
   plot3d( #lower hemisphere
      [r, theta, -sqrt(4-r^2)], r= 1..2, theta= 0..2*Pi,
      coords= cylindrical, transparency= .4
   ),
   plot3d( #inner cylinder
      [1, theta, z], theta= 0..2*Pi, z= -sqrt(3)..sqrt(3),
      coords= cylindrical, color= pink, style= patchnogrid
   ),
   scaling= constrained
);

I compute the volume five ways. In each of the triple integrals, I compute the volume in the first octant and multiply by 8.

1. As a triple integral in cylindrical coordinates:

8*Int(r, [z,r,theta]=~ [0..sqrt(4-r^2), 1..2, 0..Pi/2]);

8*(Int(r, [z = 0 .. (-r^2+4)^(1/2), r = 1 .. 2, theta = 0 .. (1/2)*Pi]))

value(%);

4*3^(1/2)*Pi

2. As a triple integral in Cartesian coordinates:

8*(Int(1, [z,y,x]=~ [0..sqrt(4-(x^2+y^2)), sqrt(1-x^2)..sqrt(4-x^2), 0..1]) +
   Int(1, [z,y,x]=~ [0..sqrt(4-(x^2+y^2)), 0..sqrt(4-x^2), 1..2])
   );

8*(Int(1, [z = 0 .. (-x^2-y^2+4)^(1/2), y = (-x^2+1)^(1/2) .. (-x^2+4)^(1/2), x = 0 .. 1]))+8*(Int(1, [z = 0 .. (-x^2-y^2+4)^(1/2), y = 0 .. (-x^2+4)^(1/2), x = 1 .. 2]))

value(%);

4*3^(1/2)*Pi

3. As a triple integral in spherical coordinates (using phi as azimuth/longitude and theta as zenith/latitude):

8*Int(rho^2*sin(theta), [rho, theta, phi]=~ [csc(theta)..2, Pi/6..Pi/2, 0..Pi/2]);

8*(Int(rho^2*sin(theta), [rho = csc(theta) .. 2, theta = (1/6)*Pi .. (1/2)*Pi, phi = 0 .. (1/2)*Pi]))

value(%);

4*3^(1/2)*Pi

4. As the volume of a solid of revolution about the y-axis by the method of cylidrical shells. (The portion in the first quadrant is rotated.)

2*Int(2*Pi*x*sqrt(4-x^2), x= 1..2);

2*(Int(2*Pi*x*(-x^2+4)^(1/2), x = 1 .. 2))

value(%);

4*3^(1/2)*Pi

5. As the volume of a solid of revolution about the x-axis by the method of washers. (The portion in the first quadrant is rotated.)

2*Int(Pi*((4-x^2)-1), x= 0..sqrt(3));

2*(Int(Pi*(-x^2+3), x = 0 .. 3^(1/2)))

value(%);

4*3^(1/2)*Pi

 


 

Download SphereCylinder.mw

You should state your Maple version in your Question header.

It works for me in Maple 2016:

restart:
gg:=A* exp( - ( (t - t0) / (tau) )**2 ):
val1:=int(gg, t=-x0..x1) assuming t0::real, tau::real, x0<x1;

val1 := (1/2)*erf((x0+t0)/tau)*A*sqrt(Pi)*tau-(1/2)*erf((-x1+t0)/tau)*A*sqrt(Pi)*tau

You should put your Maple version in your Question header.

It works for me in Maple 2016:

restart:
gg:=A* exp( - ( (t - t0) / (tau) )**2 ):
val1:=int(gg, t=-x0..x1) assuming t0::real, tau::real, x0<x1;

@John Fredsted If you get that message ("The resource that you're looking for..."), it almost certainly means that a Question has been changed to a Post (or vice versa) or the author has removed their own material and the Active Conversations index hasn't yet been updated (it takes about 15 minutes). Why it takes 15 minutes to update such a trivial index is beyond me.

People's fear of censorship on Internet fora is vastly exaggerated. I'm only aware of two cases of censorship here at MaplePrimes in my years here. One was for some explicitly racist content, and the other was for posting some copyrighted code without permission. In both cases, the offending message was replaced with a message explaining the reason for removal.

Moderators have the ability to change Posts to Questions (or vice versa) or to edit content in any other way. I often change Posts to Questions or re-attach a Reply to a new parent. Indeed, I'm tempted to change this thread to a Post.

First 199 200 201 202 203 204 205 Last Page 201 of 395