## 13202 Reputation

12 years, 350 days

## More of an observation...

because obviously the attached, using Maple's inbult numeric solvers, just *works"

So why do anything else?

However if you do (for some strange reason)  wabt to use the "differential transformation method", then I would mak the following observations

1. I know nothing about the differential transformation method
2. It appears to be related to "power series" solution methods. In fact there appears to be some controversy, in the online literature  as to whether it is actually different from a Taylor series method, or is simply a disguised variant.
3. The one thing that *all power series*  methods share is that they rely on a power series expansion about a single point. So they are appropriate for initial-value problems, where all  (so-called) boundary conditions refer to a single point - in other word they are not "boundary conditions", they are in fact "initial conditions"
4. Your problem has boundary conditions defined at eta=0 and eta=infinity - so about which of these points do you wish to perform a power series expansion?
5. I have no doubt that some here will suggest that I convert the OP's  boundary-value problem to an intial-value problem using the shooting method. Definitely possible, but I am seriously reluctant to do this when one of the boundaries is "infinity"
 > restart:   eqn1 := diff(f(eta), `\$`(eta, 3))+f(eta)*(diff(f(eta), `\$`(eta, 2)))-(diff(f(eta), eta))^2-lambda*(diff(f(eta), eta)) = 0:   eqn2 := diff(theta(eta), `\$`(eta, 2))+f(eta)*(diff(theta(eta), eta))*Pr = 0:   Bcs := (D(f))(0) = 1, f(0) = 0, (D(f))(infinity) = 0, theta(0) = 1, theta(infinity) = 0:   params:=[lambda = .5, Pr = 6.3, infinity=10]:   sol:=dsolve( eval([eqn1, eqn2, Bcs], params), numeric):   plots:-odeplot( sol,                   [ [ eta, f(eta) ],                     [ eta, theta(eta) ]                   ],                   eta=0..10,                   color=[red, blue]                 );
 >

## Still a complete...

misunderstanding of viable boundary conditions. Thje code

```with(DETools):
with(plots):
with(IntegrationTools):
Nb:=1:Nt:=1:h(x):=exp(x):
Eq1 := (diff(f(x),x,x))+(((1/x)*(diff(f(x),x))+Nb*((diff(f(x),x))*diff(g(x),x))+Nt(diff(f(x),x)^2))):
Eq2 := (diff(g(x),x,x))+(((1/x)*(diff(g(x),x))+(Nb/Nt)*((diff(f(x),x,x))+(1/x)*diff(f(x),x)))):

# Boundary condition f(x)=0, g(x)=0 at x=h(x)" where h(x) is some function of x",

#f '(x)=0, g'(x)=0 at x=0

Cd1 := f(h(x)) = 0, (D(f))(0) = 0 :
dsys := {Cd1, Eq1}:
dsol := dsolve(dsys, numeric, output = operator):
plots[odeplot](dsol, [x, diff(f(x), x\$1)], 0 .. 5, color = green):
Cd2 := g(h(x)) = 0, (D(g))(0) = 0:
dsys := {Cd1, Cd2, Eq1, Eq2}:
dsol := dsolve(dsys, numeric, output = operator):
plots[odeplot](dsol, [x, f(x)], 0 .. 5, color = red);
plots[odeplot](dsol, [x,g(x)], 0 .. 5, color = black);```

contains two *interesting* definitions, namely

```h(x):=exp(x)
f(h(x))=0```

So just consider these for a moment. For any (real) value of x from -infinity to +infinity, exp(x) varies from 0 to +infinity, so h(x) is any value from 0 to +infinity - and the you define f(h(x))=0 - so f() for any value from 0 to +infinity is zero. You think this is a valid "boundary" condition??

## As far as I know...

there is no "tutor" available for use with dsolve().

One of the stategies used by  the dsolve() command is to try to identify whether a supplied ODE is of a certain "type".

The code

```#
# Define system and solve
#
sys:=[ diff(y(x),x) = y(x)^4 * cos(x) + y(x)*tan(x),
y(0) = 0.5
]:
[_Bernoulli]

```

shows that Maple can identify the ODE as a Bernoulli equation. If you change the odeadviser() command in the above to

`DETools:-odeadvisor(sys[1], help);`

then Maple will display  help page describing the general strategy used for solving equations of this type

## With Maple...

there is very often more than one way to do something

## Rule number 1...

Use the big green up-arrow in the Mapleprimes toolbar

## You refuse...

The ode you define in your worksheets as  'de1' only contains the dependent variable 'f(x)'. The ode you define in your worksheets as  'de2' contains the dependent variables 'f(x)' and 'g(x)'.

When you process the ode 'de1', to produce 'DE1', the latter contains the dependent variables

{b[0](x), b[1](x), b[2](x), b[3](x), b[4](x), b[5](x), b[6](x), b[7](x), b[8](x), b[9](x), b[10](x)}

However when you process the ode 'de2', to produce 'DE2', the latter contains the dependent variables (check the one highlighted in red)

{f(x), c[0](x), c[1](x), c[2](x), c[3](x), c[4](x), c[5](x), c[6](x), c[7](x), c[8](x), c[9](x), c[10](x)}

The mere existence of the undetermined function 'f(x)' in the definition of 'DE2' is the source of your problem.

See the attached for proof - it makes no attempt to solve anything - it just lists the indeterminates. Until you fix the existence of the indeterminate 'f(x)' in the definition of 'DE2' you are going to get precisely nowhere

 > restart:
 > Pr:=1:
 > de1 := (1 - p)*diff(f(x), x \$ 3) + p*(diff(f(x), x \$ 3) + 1/2*f(x)*diff(f(x), x \$ 2)):
 > de2 := (1 - p)*diff(g(x), x \$ 2)/Pr + p*(diff(g(x), x \$ 2)/Pr + 1/2*f(x)*diff(g(x), x)):
 > ibvc := f(0), D(f)(0), D(f)(5) - 1, g(0) - 1, g(5):   n := 10:   F := unapply(add(b[k](x)*p^k, k = 0 .. n), x):   G := unapply(add(c[k](x)*p^k, k = 0 .. n), x):
 > DE1 := series(eval(de1, f = F), p = 0, n + 1):   indets(DE1, function(name));
 (1)
 > DE2 := series(eval(de2, g = G), p = 0, n + 1):   indets(DE2, function(name));
 (2)
 >

## Have you tried...

actually looking at the system you are trying to solve?

The final execution group in your latest worksheet - ie this piece of code

```for k from 0 to n do
IBVC1 := select(has, CT, c[k]);
{coeff(DE2, p, k), op(IBVC1)};
slv1 := dsolve({coeff(DE2, p, k), op(IBVC1)});
c[k] := unapply(rhs(slv1), x);
end do```

For k=0, your code generates the ODE system

{c[0](0) - 1, diff(c[0](x), x, x), c[0](5)}

a second order ODE with two boundary conditions - which is fine, and Maple generates the solution

slv1 := c[0](x) = -x/5 + 1

For k=1, your code generates the ODE system

{diff(c[1](x), x, x) - f(x)/10, c[1](0), c[1](5)}

a second-order ODE with two boundary conditions, and a completely unknown function of the same independent variable f(x). What solution do you expect?

Maple returns a "formal" solution containing double integrals of the unknown function f(x)

## Not sure what point...

you are trying to make.

1. Your original code generates an invalid ODE system, which means that dsolve() will produce an error.
2. This update to your code produces a valid ODE system, which means that dsolve() has a reasobale chance of coming up with a solution.

These two cases, both for loop index k=0 are illustrated in the attached. One errors, for reasons I explained previously and one doesn't.

 > # # OP's original code asks for a solution of #   dsolve({1, diff(c[0](x), x, x)}); # # which will produce an error #
 > # # OP's new code asks for a solution of #   dsolve( {D(b[0])(5) - 1, diff(b[0](x), x, x, x), b[0](0), D(b[0])(0)});
 (1)
 >

## Use the ignore=true option...

as in the attached "toy" example

 > restart;   with(Statistics):   V:=Vector[column](9, [124.0, 130.0, 130.0, 119.0, 136.0, 118.0, undefined, 130.0, 95.0]);   Mean(V, ignore=true);
 (1)
 >

## I assume...

that you were trying to achieve what is shown in the attached? I don't think this approach can be (easily) extended beyond first differences.

 > restart;   a := 1: b := 5: h := 1: f := 1/x: N := (b-a)/h:   for i from 0 while i <= N do       x[i] := h*i+a;       y[i] := eval(f, x = x[i])   end do:   printf("_______________________________________________________________________________\n\n\n");     for i from 0 by 1 while i<=2*N do:         if   irem(i,2)=0       then printf("%2d%16.7f%16.7f\n",i,x[i/2],y[i/2]);       else printf("%2d%48.7f\n",i, y[(i+1)/2]-y[(i-1)/2]);       fi;   od;   printf("_______________________________________________________________________________\n")
 _______________________________________________________________________________  0       1.0000000       1.0000000  1                                      -0.5000000  2       2.0000000       0.5000000  3                                      -0.1666667  4       3.0000000       0.3333333  5                                      -0.0833333  6       4.0000000       0.2500000  7                                      -0.0500000  8       5.0000000       0.2000000 _______________________________________________________________________________
 >

@ check my original comments on the typos in your equation and fix them properly before doing anything else.

I'd be surprised if an analytic solution can be found for the resulting ODE, so you may have to be satisfied with a numeric one - which means that you will need values for all parameters and a couple of initial/bondary conditions:-(

## With a trivial modification...

the attached will show which months in a specified year have a Friday 13th. Output is in the form

[year, [list of months with Friday 13th]]

From which is pretty obvious that in 2022, the only month with a Friday 13th is month 5 - ie May!

 > restart;   with(Calendar):   fri:= yr-> local j:seq(`if`(DayOfWeek(yr, j, 13)=6,j, NULL), j=1..12):
 > # # So which months in the supplied year have a Friday 13th #   seq( [j, [fri(j)]], j=2000..2022);
 (1)
 >
 >

## You mean...

when you type

in Maple Flow - it doesn't work?

Because if this is true then all you have demonstrated ist that Maple Flow is unable to execute basic Maple commands - which is somewhat scary!

## Idle curiosity...

If you want to restrict the range variable to integers, then not only do you need the 'sample' option, but ialso to set the adaptive=false option. If the intent is to restrict the range variable to integers, why have a 'floor(n)' command in the function definition? Isn't this a bit - well - superfluous?

## The point I was trying to make...

was that in a plot() command, the plotting variable is continuous - it does not assume integer values. Kitonum circumvented this problem by plotting only points. This has the drawback that if you draw lines between the points (say by using style=pointline), then you will not get "vertical" lines, since x-values will always differ by one

I gave a (quick+dirty) workaround for this - and you can see the effect of superimposing my plot and Kitonum's in the attached. This gives an apparent "right shift to the curve on my plot. If this is a problem it is relatively easy to fix, see the final figure in the attached

 > restart;   f:=n->ceil(sqrt(4*floor(n)))-floor(sqrt(2*floor(n)))-1:   plot(f, 10..100, gridlines=false);
 > Points:=[[10, 2], [11, 2], [12, 2], [13, 2], [14, 2], [15, 2], [16, 2], [17, 3], [18, 2], [19, 2], [20, 2], [21, 3], [22, 3], [23, 3], [24, 3], [25, 2], [26, 3], [27, 3], [28, 3], [29, 3], [30, 3], [31, 4], [32, 3], [33, 3], [34, 3], [35, 3], [36, 3], [37, 4], [38, 4], [39, 4], [40, 4], [41, 3], [42, 3], [43, 4], [44, 4], [45, 4], [46, 4], [47, 4], [48, 4], [49, 4], [50, 4], [51, 4], [52, 4], [53, 4], [54, 4], [55, 4], [56, 4], [57, 5], [58, 5], [59, 5], [60, 5], [61, 4], [62, 4], [63, 4], [64, 4], [65, 5], [66, 5], [67, 5], [68, 5], [69, 5], [70, 5], [71, 5], [72, 4], [73, 5], [74, 5], [75, 5], [76, 5], [77, 5], [78, 5], [79, 5], [80, 5], [81, 5], [82, 6], [83, 6], [84, 6], [85, 5], [86, 5], [87, 5], [88, 5], [89, 5], [90, 5], [91, 6], [92, 6], [93, 6], [94, 6], [95, 6], [96, 6], [97, 6], [98, 5], [99, 5], [100, 5]]: plots:-display( [ plot(Points, style=point, color=red, symbol=solidcircle),                   plot(f, 10..100)                 ],                 gridlines=false               );
 > g:=n->ceil(sqrt(4*floor(n+1/2)))-floor(sqrt(2*floor(n+1/2)))-1:   plots:-display( [ plot(Points, style=point, color=red, symbol=solidcircle),                     plot(g, 10..100)                   ],                   gridlines=false                 );
 >