We have a problem in solving a differential equation. We have a differential equation, say
de and an initial value condition: \(y(0)=0\) and a boundary value condition: D(y)(20)=0,
eg.
If we try to solve this equation using the dsolve command, Maple complains about the second condition. Besides, we use the numeric option of the dsolve command.
Do you know a way to solve such an equation with Maple?
Maple’s numeric dsolve is, as the on-line help describes, of use only for initial value problems. You are asking for the numerical solution of a two-point boundary value problem.
While you do not provide enough information to know exactly what type of ODE you are attempting to solve, I might have something that will be of interest.
I have recently developed a Maple implementation of the simple shooting method.
This procedure works very much like dsolve[numeric]. It must be used with
caution, as the numerical method is not always stable. A paper introducing this
procedure (shoot) and demonstrating its usage on three problems will appear in the
upcoming Special Issue of MapleTech. (These examples are nonlinear and some involve
systems.)
The on-line help for shoot, which includes two (simple) examples, is attached to the end of this message.
Anyone interested in this paper and/or procedure should send mail to: Douglas B. Meade
|\^/| Maple V Release 3 (University of South Carolina) ._|\| |/|_. Copyright (c) 1981-1994 by Waterloo Maple Software and the \ MAPLE / University of Waterloo. All rights reserved. Maple and Maple V <____ ____> are registered trademarks of Waterloo Maple Software. | Type ? for help. > read shoot; > ?shoot FUNCTION: shoot - shooting method for two-point boundary value problems CALLING SEQUENCE: shoot(deqns, vars, bc, init, options) PARAMETERS: deqns - ordinary differential equations in vars vars - variables to be solved for bc - boundary conditions init - initial values for shooting method options - optional arguments SYNOPSIS: - shoot uses the shooting method to compute approximate solutions to two-point boundary value problems for a first-order system of ordinary differential equations using the shooting method. The default output of this procedure is a list of procedures, as generated by numeric dsolve (see ?dsolve,numeric) - Each of the first four arguments must be specified as a single equation or a list or set of equations (very much like the arguments to dsolve). - The first two arguments are exactly the same as for dsolve. That is, the first argument contains the differential equations and initial conditions. The second argument contains the dependent variable(s) in the problem. - The initial conditions will include both the boundary values at one endpoint and the auxiliary initial conditions needed for the shooting method (but not the initial conditions for the sensitivity equations). The initial values for the auxiliary conditions will be given in terms of the shooting method parameters. There should be an equal number of differential equations and initial conditions. - The the third argument contains the boundary condition(s) -- as a single equation or as a list or set of equations -- at the second boundary point. These are the control equations for the iterations. - The final required argument is a single equation or list or set of equations containing the parameter values for the first iteration of the shooting method. - Any optional arguments must be in the form of keyword=value. These can be any valid optional argument to dsolve, or one of the following (with default settings listed in parentheses): 'bcerr' = bctol (Float(1,2-Digits)) 'iterr' = ittol (Float(1,2-Digits)) 'maxiter' = maxit (10) 'itmethod' = itmeth ('newton') where: - bctol is the maximum acceptable L2-error between the computed and exact boundary conditions; - ittol is the maximum acceptable L2-error between successive parameter estimates; - maxiter is the maximum number of iterations - itmeth is the root-finding method used to compute new parameter estimates (only Newton-Raphson is available) - The only known option to dsolve that is not supported is output=listprocedure. However, THERE IS NO GUARANTEE THAT ALL COMBINATIONS OF OPTIONS WILL WORK AS EXPECTED. - Each iteration of the shooting method begins by solving (by numeric dsolve) the initial value problem using the most recent set of shooting method parameters. Each boundary condition is then evaluated. If these conditions are not satisfied to the requested tolerance, one iteration of the Newton-Raphson method is used to generate a new set of shooting method parameters. - Additional information about each iteration can be obtained via the infolevel command. To see each set of parameters, use infolevel[shoot]:=1:. To include the error at the boundary, use infolevel[shoot]:=2:. Use infolevel[shoot]:=3: to see the full system of differential equations involved in the iterations. To reset to normal usage, use infolevle[shoot]:=0:. - For more information on the shooting method, see Stoer and Bulirsch (p. 469) or other standard numerical analysis texts. EXAMPLES: > restart; > with(plots): > read shoot; -------------------------------------------------------------------------------- # Example 1: x'' + x = cos(t), x(1)=3, x'(3)=1 > ODE := { diff( x(t), t ) = y(t), diff( y(t), t ) = -x(t)+cos(t) }; d d ODE := {---- x(t) = y(t), ---- y(t) = - x(t) + cos(t)} dt dt > FNS := { x(t), y(t) }; FNS := {x(t), y(t)} > IC := { x(1)=3, y(1)=alpha }; IC := {x(1) = 3, y(1) = alpha} > BC := y(3)=1; BC := y(3) = 1 > S := shoot( ODE union IC, FNS, BC, alpha=0 ); S := [t = proc(t) ... end, x(t) = proc(t) ... end, y(t) = proc(t) ... end] > odeplot( S, [t,x(t)], 1..3 ); -------------------------------------------------------------------------------- # Example 2: x''+x = cos(t), x(1)=3, x(3)+x'(3)=1 > infolevel[shoot]:=1: > S := shoot( ODE union IC, FNS, x(3)+y(3)=1, alpha=0, bcerr=Float(5,-7) ): shoot: New parameter values: alpha = 0 shoot: New parameter values: alpha = 12.08987391 shoot: New parameter values: alpha = 12.08987493 > infolevel[shoot]:=0: > odeplot( S, [ [t,x(t)], [t, x(t)+y(t)] ], 1..3, 0..20, labels=[x,t] ); # The exact solution is not hard to find for this problem: > X := rhs( dsolve( { diff(x(t),t$2)+x(t)=cos(t), x(1)=3 }, x(t), method=laplace ) ): > Xp := diff(X,t): > c := solve( subs(t=3,X+Xp)=1, {D(x)(1)} ): # Note that the final value of the shooting parameter is correct # to 5 digits to the right of the decimal point > evalf(c); {D(x)(1) = 12.08987954} > X := subs(c,X): evalf(X); .4207354924 (t - 1.) cos(t - 1.) + .2701511530 (t - 1.) sin(t - 1.) + 11.66914405 sin(t - 1.) + 3. cos(t - 1.) # The results are very comparable when viewed graphically > p1 := odeplot( S, [t,x(t)], 1..3, style=point ): > p2 := plot( X, t=1..3, color=red ): > display( {p1,p2} );
SEE ALSO: dsolve[numeric], userinfo PS: It worked very well for our problem. U.
Klein