tomleslie

13876 Reputation

20 Badges

15 years, 163 days

MaplePrimes Activity


These are answers submitted by tomleslie

because from the equations you provide, the only "sensible" conclusion I can draw is shown in the attached.

This probably isn't what you want!

But is correct, based on the information you have provided

restart

AA := [diff(ZETA[1](t), t) = -(1/3)*C[2, 1](t), diff(ZETA[1](t), t) = -(1/5)*C[4, 1](t), diff(ZETA[2](t), t) = -(1/3)*C[2, 2](t), diff(ZETA[2](t), t) = -(1/5)*C[4, 2](t), diff(ZETA[3](t), t) = -(1/3)*C[2, 3](t), diff(ZETA[3](t), t) = -(1/5)*C[4, 3](t), diff(ZETA[4](t), t) = -(1/3)*C[2, 4](t), diff(ZETA[4](t), t) = -(1/5)*C[4, 4](t)]

[diff(ZETA[1](t), t) = -(1/3)*C[2, 1](t), diff(ZETA[1](t), t) = -(1/5)*C[4, 1](t), diff(ZETA[2](t), t) = -(1/3)*C[2, 2](t), diff(ZETA[2](t), t) = -(1/5)*C[4, 2](t), diff(ZETA[3](t), t) = -(1/3)*C[2, 3](t), diff(ZETA[3](t), t) = -(1/5)*C[4, 3](t), diff(ZETA[4](t), t) = -(1/3)*C[2, 4](t), diff(ZETA[4](t), t) = -(1/5)*C[4, 4](t)]

(1)

seq(-rhs(AA[2*j-1]) = -rhs(AA[2*j]), j = 1 .. 4)

(1/3)*C[2, 1](t) = (1/5)*C[4, 1](t), (1/3)*C[2, 2](t) = (1/5)*C[4, 2](t), (1/3)*C[2, 3](t) = (1/5)*C[4, 3](t), (1/3)*C[2, 4](t) = (1/5)*C[4, 4](t)

(2)

 

Download crapSys.mw

 

   - why do you want to extract coordinates from a spacecurve?

If you can supply Maple with enough information to generate a spacecurve, then you can generate x-, y-, z-values for any parameter(s) - so why bother with the spacecurve? Sure, it will give you values, for z, y,  z, for a restricted set of the parameter(s) values which you have, where these parameter(s) values are determined by Maple's adaptive plotting routines - but why not just calculate any vakue if x, y, z. directly. Because iif the spacecurve() command can dio it, then so can you.

Unfortunately you prcvide insufficient infirmation to demonstrate how this might be done. The only things I can surmise form the completely inadequate code which you present is that each of 'x', 'y', 'z' are functions of a single variable. So for some value, say 'alpha' of this variable, why not just compute x(alpha), y(alpha), z(alpha), and forget the spacecurve altogether?

Perhaps if you had provided sufficient code to make your question "executable", at least capable of producing the spacecurve, then I could make this explanation clearer.

Rule No 1: Ask an unclear question and you will get an unclear answer

minor typos, which are fixed in the attached

  restart:
  with(plots):
  with(plottools):
  Vdot := proc (U, V)local i: add(U[i]*V[i], i = 1 .. 2) end proc:
  dist := proc (M, N) sqrt(Vdot(expand(M-N), expand(M-N))) end proc:
  ngon := n -> local i: [seq([cos(2*Pi*i/n), sin(2*Pi*i/n)], i = 1 .. n)]:
  theta := (2*Pi)/5:
  poly := [seq([cos(k*theta), sin(k*theta)], k = 1 .. 5)]:
  Ii := [0, 1/2]:
  H := [-1/4, 0]:
  r := dist(Ii, H):
  theta := (2*Pi)/5:
  p1 := pointplot([seq([cos(k*theta), sin(k*theta)], k = 0 .. 5)], symbol = solidcircle, color = red, symbolsize = 10):
  p2 := textplot([seq([cos(k*theta), sin(k*theta), cat("M", k)], k = 0 .. 4)], align = ["above", "right"]):
  cir1 := circle([0, 0], 1/2, color = green, linestyle=dashdot):
  cir2 := circle([-1/4, 0], r, color = black):
  cir3:=circle([0,0],1,color=red):
  display([p1, p2, cir1, cir2, cir3, polygonplot(poly, thickness = 5, color = blue, transparency = 0.95)], axes = normal);

 

 

Download pentplot.mw

It appears that visualization routines within "Student" packages come with their own set of plot options, whch are only a subset of those available for "normal" plot() and plot3d() commands. You can see the available options using ?Student Plot Options.

To my amazement, no "color" option is available!!!. So in the OP's case, in the option,

volumeoptions = [color = blue, transparency = 0.6]

The color sub-option is simply ignored !

It seems that the visualization routines within student packages come with a default palette - so the first object in a plot will always use the first color in the palette, second object will be the second palette color, and so on. By default, the color palette is "Niagara" which you can view here ?Niagara Color Palette.

Thus the only way to change colors is to add desired colors to the beginning of the default palette, using the Student:-SetColors() command, whihch I have done in the attached to add blue, red green. Note that in the resulting plot, the first object (ie the volume) is "blue" and the second object (ie the generating line) is red.

  restart:
  Student:-SetColors(blue, red, green);
  with( Student:-Calculus1):
  VolumeOfRevolution( 4,
                      x = 1 .. 2,
                      volumeoptions = [transparency = 0.6],
                      orientation = [270, 0, 15],
                      view = [0 .. 3, -5 .. 5, -5 .. 5],
                      labels = [x, y, z],
                      output = plot,
                      axis = horizontal,
                      scaling = constrained
                    );

[blue, red, green, "#3E578A", "#780072", "#00786A", "#604191", "#004A78", "#784C00", "#91414A", "#3E738A", "#78003B", "#00783F", "#914186", "#510078", "#777800"]

 

 

 

Download addColStud.mw

 

They are two entirely different concepts.

If you need a "toy" example to illustrate the difference, consider sin(x)/x. This has a limiting value of 1 as x->0. Hoever if you substitute x=0, you get 0/0 - which is what, exactly?

This site really isn't the best place to learn high-school mathematics

you state

Here psi is a general wave function from schrodinger wave equation.

In which case it probably ought to be Psi(x,t) - but I wouldn't use that name, because Psi() is one of Maple's built-in functions and Psi(x,t) will be interpreted as  the x-th derivative of the digamma function Psi(t). Just avoid the problem by using Phi(x,t) instead,

In the attached I have unapply() to evaluate the arguments of all necessary functions as far as possible.

One still doesn't get too far, because there are still many unknown functions in the final integrand eg Phi(x,t), mu(t), and x(t) - so no integration can actually be performed.

Another piece of advice I wouldn't use a functon 'x(t)' and a scalar independent variable named 'x' in the same worksheet - this is just asking for "subtle" problems.

Anyhow, for what it is worth, you may(?) find the attached helpful

  restart:

  f := unapply
       ( abs(Phi(x,t))^2,
         [x,t]
       );

proc (x, t) options operator, arrow; abs(Phi(x, t))^2 end proc

(1)

  a := unapply
       ( piecewise
         ( 0 <= t and t <= 1,
           1.5*t,
           1 <= t and t <= 2,
           1.5*(2 - t)
         ),
         t
       );

a := proc (t) options operator, arrow; piecewise(0 <= t and t <= 1, 1.5*t, 1 <= t and t <= 2, 3.0+(-1)*1.5*t) end proc

(2)

  y[a] := unapply
          ( piecewise
            ( 0 <= t and t <= 0.1,
              a(t),
              0.1 <= t and t <= 0.2,
              -a(t)
            ),
            t
          );

y[a] := proc (t) options operator, arrow; piecewise(0 <= t and t <= .1, piecewise(0 <= t and t <= 1, 1.5*t, 1 <= t and t <= 2, 3.0+(-1)*1.5*t), .1 <= t and t <= .2, -piecewise(0 <= t and t <= 1, 1.5*t, 1 <= t and t <= 2, 3.0+(-1)*1.5*t)) end proc

(3)

  y := unapply
       ( y[a](t) + mu(t),
         t
       );

y := proc (t) options operator, arrow; piecewise(0 <= t and t <= .1, piecewise(0 <= t and t <= 1, 1.5*t, 1 <= t and t <= 2, 3.0+(-1)*1.5*t), .1 <= t and t <= .2, -piecewise(0 <= t and t <= 1, 1.5*t, 1 <= t and t <= 2, 3.0+(-1)*1.5*t))+mu(t) end proc

(4)

  w := unapply
       ( int
         ( x(t)*f(x, t),
           x
         ),
         t
       );

proc (t) options operator, arrow; int(x(t)*abs(Phi(x, t))^2, x) end proc

(5)

  v := unapply
       ( y(t) - w(t)*w(t),
         t
       );

v := proc (t) options operator, arrow; piecewise(0 <= t and t <= .1, piecewise(0 <= t and t <= 1, 1.5*t, 1 <= t and t <= 2, 3.0+(-1)*1.5*t), .1 <= t and t <= .2, -piecewise(0 <= t and t <= 1, 1.5*t, 1 <= t and t <= 2, 3.0+(-1)*1.5*t))+mu(t)-(int(x(t)*abs(Phi(x, t))^2, x))^2 end proc

(6)

  eqn:= diff(K(x, t), t) = beta*v(t)*f(x, t);

diff(K(x, t), t) = beta*(piecewise(0 <= t and t <= .1, piecewise(0 <= t and t <= 1, 1.5*t, 1 <= t and t <= 2, 3.0-1.5*t), .1 <= t and t <= .2, -piecewise(0 <= t and t <= 1, 1.5*t, 1 <= t and t <= 2, 3.0-1.5*t))+mu(t)-(int(x(t)*abs(Phi(x, t))^2, x))^2)*abs(Phi(x, t))^2

(7)

  sol:= map( int, eqn, t);

K(x, t) = int(beta*(piecewise(0 <= t and t <= .1, piecewise(0 <= t and t <= 1, 1.5*t, 1 <= t and t <= 2, 3.0-1.5*t), .1 <= t and t <= .2, -piecewise(0 <= t and t <= 1, 1.5*t, 1 <= t and t <= 2, 3.0-1.5*t))+mu(t)-(int(x(t)*abs(Phi(x, t))^2, x))^2)*abs(Phi(x, t))^2, t)

(8)

#
# Expand the integral on the rhs of the above
#

  lhs(sol)=IntegrationTools:-Expand(rhs(sol));

K(x, t) = beta*(int(abs(Phi(x, t))^2*piecewise(0 <= t and t <= .1, piecewise(0 <= t and t <= 1, 1.5*t, 1 <= t and t <= 2, 3.0-1.5*t), .1 <= t and t <= .2, -piecewise(0 <= t and t <= 1, 1.5*t, 1 <= t and t <= 2, 3.0-1.5*t)), t))+beta*(int(abs(Phi(x, t))^2*mu(t), t))-beta*(int(abs(Phi(x, t))^2*(int(x(t)*abs(Phi(x, t))^2, x))^2, t))

(9)

 

Download oddint.mw

 

I think this might depend (a lot?) on how "generally applicable" you want the process to be. For the example you give,the simple code in the attached will work, For other (possibky similar, but slightly different?) problems, some "tuning" might be necessary, but I wouldn't expect such "tuning" to be difficult

  with(numapprox):
  op~(1, [op(chebyshev(exp(x),x))]);

[HFloat(1.2660658777520082), HFloat(1.1303182079849698), HFloat(0.2714953395340765), HFloat(0.044336849848663804), HFloat(0.005474240442093705), HFloat(5.429263119139931e-4), HFloat(4.497732295427603e-5), HFloat(3.1984364624457973e-6), HFloat(1.992124806415823e-7), HFloat(1.103677170949997e-8), HFloat(5.505896979790788e-10)]

(1)

 

Download extrCoeff.mw

produces correct results! You can see this by reconstructing the original equations by using the matrix and vector of right-hand-sides produced by the GenerateMatrix() command to reconstruct the original equations. This is shown in the final execution group of the attached.


(Normally one would do this "reconstruction" using the GenerateEquations() command but this doesn't seems to like the fact that your designated "variables" are not simple names.)
 

restart

N := 2; M := 2

a := 1b := 1

Qa := [-0.5379667864e-1*(diff(tau[1, 1](t), t, t))+7.862351349*10^(-11)*tau[2, 1](t)-8.050993899*10^(-12)*(diff(tau[2, 1](t), t, t))+.1166068042*(diff(tau[1, 2](t), t))+2.181309895*10^(-11)*(diff(tau[2, 2](t), t))+.5309519363*tau[1, 1](t) = 0, -1.265965258*10^(-11)*(diff(tau[1, 1](t), t, t))+.4884414390*tau[2, 1](t)-0.4948946475e-1*(diff(tau[2, 1](t), t, t))+2.738892495*10^(-11)*(diff(tau[1, 2](t), t))+.1340883970*(diff(tau[2, 2](t), t))+1.246469610*10^(-10)*tau[1, 1](t) = 0, 3.649366137*10^(-10)*tau[2, 2](t)-9.135908950*10^(-12)*(diff(tau[2, 2](t), t, t))-5.160677740*10^(-11)*(diff(tau[2, 1](t), t))+1.953765755*tau[1, 2](t)-0.4948946473e-1*(diff(tau[1, 2](t), t, t))-.3476543209*(diff(tau[1, 1](t), t)) = 0, 2.246672656*tau[2, 2](t)-0.5690888318e-1*(diff(tau[2, 2](t), t, t))-.3198194887*(diff(tau[2, 1](t), t))+4.602903411*10^(-10)*tau[1, 2](t)-1.159417294*10^(-11)*(diff(tau[1, 2](t), t, t))-8.175817372*10^(-11)*(diff(tau[1, 1](t), t)) = 0]

Q1 := [seq(seq(diff(tau[i, j](t), t), i = 1 .. M), j = 1 .. N)]

with(LinearAlgebra)

C, R := GenerateMatrix(simplify(Qa), Q1); `~`[`=`](`~`[`-`](convert(C.Matrix(4, 1, Q1), list), convert(R, list)), 0); is(% = Qa)

[-0.5379667864e-1*(diff(diff(tau[1, 1](t), t), t))+0.7862351349e-10*tau[2, 1](t)-0.8050993899e-11*(diff(diff(tau[2, 1](t), t), t))+.1166068042*(diff(tau[1, 2](t), t))+0.2181309895e-10*(diff(tau[2, 2](t), t))+.5309519363*tau[1, 1](t) = 0, -0.1265965258e-10*(diff(diff(tau[1, 1](t), t), t))+.4884414390*tau[2, 1](t)-0.4948946475e-1*(diff(diff(tau[2, 1](t), t), t))+0.2738892495e-10*(diff(tau[1, 2](t), t))+.1340883970*(diff(tau[2, 2](t), t))+0.1246469610e-9*tau[1, 1](t) = 0, 0.3649366137e-9*tau[2, 2](t)-0.9135908950e-11*(diff(diff(tau[2, 2](t), t), t))-0.5160677740e-10*(diff(tau[2, 1](t), t))+1.953765755*tau[1, 2](t)-0.4948946473e-1*(diff(diff(tau[1, 2](t), t), t))-.3476543209*(diff(tau[1, 1](t), t)) = 0, 2.246672656*tau[2, 2](t)-0.5690888318e-1*(diff(diff(tau[2, 2](t), t), t))-.3198194887*(diff(tau[2, 1](t), t))+0.4602903411e-9*tau[1, 2](t)-0.1159417294e-10*(diff(diff(tau[1, 2](t), t), t))-0.8175817372e-10*(diff(tau[1, 1](t), t)) = 0]

 

true

(1)

 


 

Download genMat1.mw

as in the attached

  restart;
  with(plots):
  P:= [.6286420119, -.6286420119, 0., 0., 0., 0., 0., 0., 0., 0., 0.]:
  Q:= [2.106333379, 2.106333379, 4.654463885, 7.843624703, 10.99193295, 14.13546782, 17.27782732, 20.41978346, 23.56157073, 26.70327712, 29.84494078]:
  W:= zip((x,y)->[x,y], P, Q);
  pointplot( W,
             axes=boxed,
             view=[-1..1, 0..30],
             symbol=solidcircle,
             symbolsize=16,
             color=red);
  display( seq(pointplot(W[j]), j = 1 .. numelems(W)),
           insequence = true,
           axes=boxed,
           view=[-1..1, 0..30],
           symbol=solidcircle,
           symbolsize=16,
           color=red
         );
  display( seq(pointplot(W[1..j]), j = 1 .. numelems(W)),
           insequence = true,
           axes=boxed,
           view=[-1..1, 0..30],
           symbol=solidcircle,
           symbolsize=16,
           color=red
         );

[[.6286420119, 2.106333379], [-.6286420119, 2.106333379], [0., 4.654463885], [0., 7.843624703], [0., 10.99193295], [0., 14.13546782], [0., 17.27782732], [0., 20.41978346], [0., 23.56157073], [0., 26.70327712], [0., 29.84494078]]

 

 

 

 

 

 


 

Download anims.mw

If you read the help at ?pdsolve/numeric, you will find the following interesting paragraphs

The pdsolve/numeric routine uses finite difference methods to obtain these numerical solutions.

 

The solver has two modes of operation.
  - The first mode of operation uses the default method, which is a centered implicit scheme, and is capable of finding solutions for higher order PDE or PDE systems. The PDE systems must be sufficiently close to a standard form for the method to find the numerical solution.
- The second mode of operation is strictly educational, and allows specification of a particular method to solve a PDE. This mode is restricted to a single PDE.

So the built-in command is just pdsolve( pdsystem, numeric) - which will use a centered implicit finite difference method.

I assume that codes from the textbook you specified are avaliable from its associated website, which appears to be here

https://sites.google.com/site/numericalanalysis1burden

  1. When you want to use the trigonometric functions sin() and cos(), the syntax is sin(alpha) and cos(alpha). The parentheses are mandatory
  2. You have multiple instances of "missing" multiplication signs. If you insist on using 2-D input mode, then you must first understand the conditions under which (so-called) implied multiplication will be assumed - and when it won't. In the attached, I fixed these by guesswork, but I wouldn't bet that I got your intention correct in all cases
  3. Don't use indexed variables such as K[f] and rho[s], when all you want is an inert subscript: use K__f and rho__s instead.
  4. When using PDETools:-declare(), 'f(x)' will display as 'f'. This does not mean" that you can use 'f' on entry and expect it to be interprted as f(x). You do this in eq2, which contains the quantity 'f'. Unless this is intended to be f(x), your equations are completely uncoupled and can be solved independently (provided you give a value for the parameter 'f', which is completely unrelated to the function f(x)
  5. You use two parameters 'phi' and 'alpha', whihc are nowhere assigned values - so in the attached I just guessed and set alpha=Pi/4, phi=0.1. You will need to change these.
  6. Although you have a variable named 'bcs', this is a siimple sequence of parameter values, it contains no boundary and/or initial conditions. Your ODEs are third -order in f(x) and second-order in Theta(x). So for numeric solution you need three boubndary conditions on f(x) and two on Theta(x). In the attached, I used f(0)=1, D(f)(0)=0,(D@@2)(f)(0)=0, Theta(0)=0, D(Theta)(0)=1, which of course are just random guesses.

With all of the above "fixes", some of which are certain to be incorrect the attached will execute and produce an answer.

restart:
with(PDETools, declare):
with(plots):
declare(f(x));

f(x)*`will now be displayed as`*f

(1)

eq1 := (diff(f(x), `$`(x, 3)))/((1-phi)^2.5*(1-phi+phi*rho__s/rho__f))+(1/2)*cos(alpha)*eta*(diff(f(x), `$`(x, 2)))+(1/2)*sin(alpha)*f(x)*(diff(f(x), `$`(x, 2)))-M*(diff(f(x), `$`(x, 1))-1)/(1-phi+phi*rho__s/rho__f) = 0;

(diff(diff(diff(f(x), x), x), x))/((1-phi)^2.5*(1-phi+phi*rho__s/rho__f))+(1/2)*cos(alpha)*eta*(diff(diff(f(x), x), x))+(1/2)*sin(alpha)*f(x)*(diff(diff(f(x), x), x))-M*(diff(f(x), x)-1)/(1-phi+phi*rho__s/rho__f) = 0

(2)

eq2 := K__nf*(diff(Theta(x), `$`(x, 2)))/K__f+(1/2)*Pr*(1-phi+phi*rho__s*c__s/(rho__f*c__f))*(eta*cos(alpha)+f(x)*sin(alpha))*(diff(Theta(x), `$`(x, 1))) = 0;

 

 

K__nf*(diff(diff(Theta(x), x), x))/K__f+(1/2)*Pr*(1-phi+phi*rho__s*c__s/(rho__f*c__f))*(eta*cos(alpha)+f(x)*sin(alpha))*(diff(Theta(x), x)) = 0

(3)

bcs := eta = 1, rho__s = 5200, rho__f = 997.1, c__s = 670, c__f = 4179, K__s = 6, K__f = .613, K__nf = K__f*(K__s+2*K__f-2*phi*(K__f-K__s))/(K__s+2*K__f+phi*(K__f-K__s));

eta = 1, rho__s = 5200, rho__f = 997.1, c__s = 670, c__f = 4179, K__s = 6, K__f = .613, K__nf = K__f*(K__s+2*K__f-2*phi*(K__f-K__s))/(K__s+2*K__f+phi*(K__f-K__s))

(4)

a1 := [cos(beta) = 1, sin(beta) = 0, M = 1, Pr = 6.2, phi = .1];
a2 := [cos(beta) = 0, sin(beta) = 1, M = 1, Pr = 6.2, phi = .1];
a3 := [cos(beta) = 1/2, sin(beta) = sqrt(3)/2, M = 1, Pr = 6.2, phi = .1];
a4 := [cos(beta) = 1, sin(beta) = 0, M = 1, Pr = 6.2, phi = .1];
a5 := [cos(beta) = 0, sin(beta) = 1, M = 1, Pr = 6.2, phi = .1];
a6 := [cos(beta) = 1/2, sin(beta) = sqrt(3)/2, M = 1, Pr = 6.2, phi = .1];

[cos(beta) = 1, sin(beta) = 0, M = 1, Pr = 6.2, phi = .1]

 

[cos(beta) = 0, sin(beta) = 1, M = 1, Pr = 6.2, phi = .1]

 

[cos(beta) = 1/2, sin(beta) = (1/2)*3^(1/2), M = 1, Pr = 6.2, phi = .1]

 

[cos(beta) = 1, sin(beta) = 0, M = 1, Pr = 6.2, phi = .1]

 

[cos(beta) = 0, sin(beta) = 1, M = 1, Pr = 6.2, phi = .1]

 

[cos(beta) = 1/2, sin(beta) = (1/2)*3^(1/2), M = 1, Pr = 6.2, phi = .1]

(5)

b1 := subs(a1, eq1);
b2 := subs(a2, eq1);
b3 := subs(a3, eq1);
b4 := subs(a4, eq1);
b5 := subs(a5, eq1);
b6 := subs(a6, eq1);

1.301348831*(diff(diff(diff(f(x), x), x), x))/(.9+.1*rho__s/rho__f)+(1/2)*cos(alpha)*eta*(diff(diff(f(x), x), x))+(1/2)*sin(alpha)*f(x)*(diff(diff(f(x), x), x))-(diff(f(x), x)-1)/(.9+.1*rho__s/rho__f) = 0

 

1.301348831*(diff(diff(diff(f(x), x), x), x))/(.9+.1*rho__s/rho__f)+(1/2)*cos(alpha)*eta*(diff(diff(f(x), x), x))+(1/2)*sin(alpha)*f(x)*(diff(diff(f(x), x), x))-(diff(f(x), x)-1)/(.9+.1*rho__s/rho__f) = 0

 

1.301348831*(diff(diff(diff(f(x), x), x), x))/(.9+.1*rho__s/rho__f)+(1/2)*cos(alpha)*eta*(diff(diff(f(x), x), x))+(1/2)*sin(alpha)*f(x)*(diff(diff(f(x), x), x))-(diff(f(x), x)-1)/(.9+.1*rho__s/rho__f) = 0

 

1.301348831*(diff(diff(diff(f(x), x), x), x))/(.9+.1*rho__s/rho__f)+(1/2)*cos(alpha)*eta*(diff(diff(f(x), x), x))+(1/2)*sin(alpha)*f(x)*(diff(diff(f(x), x), x))-(diff(f(x), x)-1)/(.9+.1*rho__s/rho__f) = 0

 

1.301348831*(diff(diff(diff(f(x), x), x), x))/(.9+.1*rho__s/rho__f)+(1/2)*cos(alpha)*eta*(diff(diff(f(x), x), x))+(1/2)*sin(alpha)*f(x)*(diff(diff(f(x), x), x))-(diff(f(x), x)-1)/(.9+.1*rho__s/rho__f) = 0

 

1.301348831*(diff(diff(diff(f(x), x), x), x))/(.9+.1*rho__s/rho__f)+(1/2)*cos(alpha)*eta*(diff(diff(f(x), x), x))+(1/2)*sin(alpha)*f(x)*(diff(diff(f(x), x), x))-(diff(f(x), x)-1)/(.9+.1*rho__s/rho__f) = 0

(6)

c1 := subs(a1, eq2);
c2 := subs(a2, eq2);
c3 := subs(a3, eq2);
c4 := subs(a4, eq2);
c5 := subs(a5, eq2);
c6 := subs(a6, eq2);

K__nf*(diff(diff(Theta(x), x), x))/K__f+3.100000000*(.9+.1*rho__s*c__s/(rho__f*c__f))*(eta*cos(alpha)+f(x)*sin(alpha))*(diff(Theta(x), x)) = 0

 

K__nf*(diff(diff(Theta(x), x), x))/K__f+3.100000000*(.9+.1*rho__s*c__s/(rho__f*c__f))*(eta*cos(alpha)+f(x)*sin(alpha))*(diff(Theta(x), x)) = 0

 

K__nf*(diff(diff(Theta(x), x), x))/K__f+3.100000000*(.9+.1*rho__s*c__s/(rho__f*c__f))*(eta*cos(alpha)+f(x)*sin(alpha))*(diff(Theta(x), x)) = 0

 

K__nf*(diff(diff(Theta(x), x), x))/K__f+3.100000000*(.9+.1*rho__s*c__s/(rho__f*c__f))*(eta*cos(alpha)+f(x)*sin(alpha))*(diff(Theta(x), x)) = 0

 

K__nf*(diff(diff(Theta(x), x), x))/K__f+3.100000000*(.9+.1*rho__s*c__s/(rho__f*c__f))*(eta*cos(alpha)+f(x)*sin(alpha))*(diff(Theta(x), x)) = 0

 

K__nf*(diff(diff(Theta(x), x), x))/K__f+3.100000000*(.9+.1*rho__s*c__s/(rho__f*c__f))*(eta*cos(alpha)+f(x)*sin(alpha))*(diff(Theta(x), x)) = 0

(7)

  sys:=eval[recurse]([b1,c1, f(0)=1, D(f)(0)=0,(D@@2)(f)(0)=0, Theta(0)=0, D(Theta)(0)=1], [bcs, alpha=Pi/4, phi=0.1 ]);
  dsol:=dsolve(sys, numeric):
  odeplot( dsol, [[ x, f(x)], [x, Theta(x)]], x=0..3, color=[red, blue]);

[.9154678101*(diff(diff(diff(f(x), x), x), x))+(1/4)*2^(1/2)*(diff(diff(f(x), x), x))+(1/4)*2^(1/2)*f(x)*(diff(diff(f(x), x), x))-.7034761075*(diff(f(x), x))+.7034761075 = 0, 1.241667040*(diff(diff(Theta(x), x), x))+3.049196273*((1/2)*2^(1/2)+(1/2)*f(x)*2^(1/2))*(diff(Theta(x), x)) = 0, f(0) = 1, (D(f))(0) = 0, ((D@@2)(f))(0) = 0, Theta(0) = 0, (D(Theta))(0) = 1]

 

 

  

 

Download badODEsys.mw

because it is difficult (for me at least) to work out exactly which form of output you want.

A couple of options are shown iin the attached.

#
# One possible variant of OP's desired expression
#
  Set:= `%+`(1/3%*Vector([3^(1/2), -3^(1/2), 3^(1/2), 0]), -1/3%*(Vector([3^(1/2), 3^(1/2), 3^(1/2), 0])+3*Vector([2,1,1,0])));
#
# Quick check on export+import to/from MathML
#
  a:= MathML:-Export(Set):
  MathML:-Import(a);

`%+`(`%*`(1/3, Vector(4, {(1) = sqrt(3), (2) = -sqrt(3), (3) = sqrt(3), (4) = 0})), -`%*`(1/3, Vector(4, {(1) = sqrt(3)+6, (2) = sqrt(3)+3, (3) = sqrt(3)+3, (4) = 0})))

 

`%+`(`%*`(1/3, Vector[column](%id = 36893488148299434812)), -`%*`(1/3, Vector[column](%id = 36893488148299434932)))

(1)

#
# Another possible variant of OP's desired expression
#
  Set:= `%+`(1/3%*Vector([3^(1/2), -3^(1/2), 3^(1/2), 0]), -1/3%*Vector([3^(1/2), 3^(1/2), 3^(1/2), 0]),Vector([2,1,1,0]));
#
# Quick check on  export+import to/from MathML
#
  a:= MathML:-Export(Set):
  MathML:-Import(a);

`%+`(`%*`(1/3, Vector(4, {(1) = sqrt(3), (2) = -sqrt(3), (3) = sqrt(3), (4) = 0})), -`%*`(1/3, Vector(4, {(1) = sqrt(3), (2) = sqrt(3), (3) = sqrt(3), (4) = 0})), Vector(4, {(1) = 2, (2) = 1, (3) = 1, (4) = 0}))

 

`%+`(`%*`(1/3, Vector[column](%id = 36893488148299426260)), -`%*`(1/3, Vector[column](%id = 36893488148299426380)), Vector[column](%id = 36893488148299426500))

(2)

 


 

Download inert2.mw

is shown in the atttached.

It seems like you are trying to do it in some complicated fashion which I don't understand.

MAybe you can tell me what you don't like abouot the attached?

  restart:
  with(Fractals:-LSystem):
  with(LSystemExamples):
  with(plots):
  display([seq( PlotExample(DragonCurve,j),j=1..12)], insequence=true);
  display([seq( PlotExample(GosperCurve,j),j=1..4)], insequence=true);

 

 

 

Download fracAnim.mw

 

why you would ever want to do this!!

But if you just use the prefix forms of the basic mathematical operations, then it is not difficult to achieve. See the attached

  restart;
  A:= Matrix(3,3, [[1,1,1],
                   [0,1,1],
                   [0,0,1]
                  ]
            ):
  F:= proc(oper)
           local J:=0, n:=3, i, j;
           for i to n do
               for j from i + 1 to n do
                   if   A[i, j] = 1
                   then J := J + abs( oper(10,20) );
                   fi;
               od;
           od;
           return J;
      end proc:
  F(`*`);
  F(`+`);
  F(`-`);
  F(`/`);

600

 

90

 

30

 

3/2

(1)

 

Download ops.mw

What is wrong with something as simple as the attached

  restart:
  params:=[ a=1, b=2, c=1, d=6, e=2]:
  P:= a*b*c*d*e;
  eval(P, params);

a*b*c*d*e

 

24

(1)

 

Download eval.mw

First 22 23 24 25 26 27 28 Last Page 24 of 207