tomleslie

13776 Reputation

19 Badges

14 years, 64 days

MaplePrimes Activity


These are answers submitted by tomleslie

I'd read the whole data set from EXxcel using ExcelTools:-Import() and then use Maple to select/organise dependent and independent variables for fitting.

A "toy" example is shown in the attached, where each fit outputs [rsquared, rsquaredadjusted, tprobability]. from the help page

?Solutions Returned by the Regression Commands

one has

  • rsquared -- The coefficient of determination, often used to indicate how well data fits a model.
  • rsquaredadjusted -- The adjusted coefficient of determination that accounts for the number of variables in the model.
  • tprobability -- The p-value for the hypothesis test for which the t value is the test statistic.

This ought(?) to be the answer to the second question which you have asked.

 

In this "toy" example the values for all of the dependent and independent variables are fictitious, as is the model to which they are being fitted so don't expect fits to be much good

By the way if you really have 30 columns of independent variables and you want to perform a fit to any pair of columns  - this will mean performing 870 fit operations -you may want to deal with a "smaller" example whilst writing/debugging your code

  restart:
  with(LinearAlgebra):
  ExperimentalData := < v, 0.531, 0.341, 0.163, 0.641, 0.713, -0.040
                      | x, 1, 1, 1, 2, 2, 2
                      | y, 1, 2, 3, 1, 2, 3
                      | z, 1, 2, 3, 4, 5, 6
                      >;
  for j in combinat:-permute(op([1,2], ExperimentalData)-1, 2) do
      labs:=`<|>`(Column(ExperimentalData[1,..], [(j+~1)[], 1]));
      data:=`<|>`(Column(ExperimentalData[2..-1,..], [(j+~1)[], 1]));
      Statistics:-LinearFit( labs[1,1]*a + b*labs[1,1]^2/labs[1,2],
                             data,
                             convert(labs[1, 1..2], list),
                             output=[rsquared, rsquaredadjusted, tprobability]
                           );
  od;

Matrix(7, 4, {(1, 1) = v, (1, 2) = x, (1, 3) = y, (1, 4) = z, (2, 1) = .531, (2, 2) = 1, (2, 3) = 1, (2, 4) = 1, (3, 1) = .341, (3, 2) = 1, (3, 3) = 2, (3, 4) = 2, (4, 1) = .163, (4, 2) = 1, (4, 3) = 3, (4, 4) = 3, (5, 1) = .641, (5, 2) = 2, (5, 3) = 1, (5, 4) = 4, (6, 1) = .713, (6, 2) = 2, (6, 3) = 2, (6, 4) = 5, (7, 1) = -0.40e-1, (7, 2) = 2, (7, 3) = 3, (7, 4) = 6})

 

Vector[row](3, {(1) = x, (2) = y, (3) = v})

 

Matrix(6, 3, {(1, 1) = 1, (1, 2) = 1, (1, 3) = .531, (2, 1) = 1, (2, 2) = 2, (2, 3) = .341, (3, 1) = 1, (3, 2) = 3, (3, 3) = .163, (4, 1) = 2, (4, 2) = 1, (4, 3) = .641, (5, 1) = 2, (5, 2) = 2, (5, 3) = .713, (6, 1) = 2, (6, 2) = 3, (6, 3) = -0.40e-1})

 

[.730248825887470, .595373238831206, [.590040726626360, .380105705667960]]

 

Vector[row](3, {(1) = x, (2) = z, (3) = v})

 

Matrix(6, 3, {(1, 1) = 1, (1, 2) = 1, (1, 3) = .531, (2, 1) = 1, (2, 2) = 2, (2, 3) = .341, (3, 1) = 1, (3, 2) = 3, (3, 3) = .163, (4, 1) = 2, (4, 2) = 4, (4, 3) = .641, (5, 1) = 2, (5, 2) = 5, (5, 3) = .713, (6, 1) = 2, (6, 2) = 6, (6, 3) = -0.40e-1})

 

[.821717463177158, .732576194765737, [.668002244228826, .133784018957921]]

 

Vector[row](3, {(1) = y, (2) = x, (3) = v})

 

Matrix(6, 3, {(1, 1) = 1, (1, 2) = 1, (1, 3) = .531, (2, 1) = 2, (2, 2) = 1, (2, 3) = .341, (3, 1) = 3, (3, 2) = 1, (3, 3) = .163, (4, 1) = 1, (4, 2) = 2, (4, 3) = .641, (5, 1) = 2, (5, 2) = 2, (5, 3) = .713, (6, 1) = 3, (6, 2) = 2, (6, 3) = -0.40e-1})

 

[.533423416834850, .300135125252275, [.142994934540483, .281965124452022]]

 

Vector[row](3, {(1) = y, (2) = z, (3) = v})

 

Matrix(6, 3, {(1, 1) = 1, (1, 2) = 1, (1, 3) = .531, (2, 1) = 2, (2, 2) = 2, (2, 3) = .341, (3, 1) = 3, (3, 2) = 3, (3, 3) = .163, (4, 1) = 1, (4, 2) = 4, (4, 3) = .641, (5, 1) = 2, (5, 2) = 5, (5, 3) = .713, (6, 1) = 3, (6, 2) = 6, (6, 3) = -0.40e-1})

 

[.380215523936068, 0.703232859041014e-1, [.405926630635241, .698700881346495]]

 

Vector[row](3, {(1) = z, (2) = x, (3) = v})

 

Matrix(6, 3, {(1, 1) = 1, (1, 2) = 1, (1, 3) = .531, (2, 1) = 2, (2, 2) = 1, (2, 3) = .341, (3, 1) = 3, (3, 2) = 1, (3, 3) = .163, (4, 1) = 4, (4, 2) = 2, (4, 3) = .641, (5, 1) = 5, (5, 2) = 2, (5, 3) = .713, (6, 1) = 6, (6, 2) = 2, (6, 3) = -0.40e-1})

 

[.917247152618675, .875870728928013, [0.583110431704181e-2, 0.968590476544140e-2]]

 

Vector[row](3, {(1) = z, (2) = y, (3) = v})

 

Matrix(%id = 36893488148279964604)

 

[HFloat(0.541823133651818), HFloat(0.312734700477727), [HFloat(0.995893258860949), HFloat(0.4743896868797166)]]

(1)

 

Download fitdata.mw

where  some typos are fixed.

I don't think this code will *ever* produce a meaningful answer, because there is no constraint on x[4], so it can be increased arbitrarily, with corresponding increases in the values of f() and g()

restart

f := proc (x) options operator, arrow; sum(ln(x[i]), i = 1 .. 4) end proc

proc (x) options operator, arrow; sum(ln(x[i]), i = 1 .. 4) end proc

(1)

g := proc (x) options operator, arrow; ln(x[1])+ln(x[2])+ln(x[3])+ln(x[4]) end proc

proc (x) options operator, arrow; ln(x[1])+ln(x[2])+ln(x[3])+ln(x[4]) end proc

(2)

f([1, 1, 1, 1])

0

(3)

cons := {x[1]+x[2] <= 1, x[1]+x[3] <= 2}

{x[1]+x[2] <= 1, x[1]+x[3] <= 2}

(4)

with(Optimization); ip := [x[1] = .5, x[2] = .5, x[3] = 1, x[4] = 100]

[x[1] = .5, x[2] = .5, x[3] = 1, x[4] = 100]

(5)

Maximize(f(x), cons, assume = nonnegative, initialpoint = ip)

[16.4237121534987303, [x[1] = HFloat(0.5028668384263846), x[2] = HFloat(0.49713316157361515), x[3] = HFloat(1.4971331615736156), x[4] = HFloat(3.62694858829851e7)]]

(6)

Maximize(f(x), {x[1]+x[2] <= 1, x[1]+x[3] <= 2}, assume = nonnegative, initialpoint = ip)

[16.4237121534987303, [x[1] = HFloat(0.5028668384263846), x[2] = HFloat(0.49713316157361515), x[3] = HFloat(1.4971331615736156), x[4] = HFloat(3.62694858829851e7)]]

(7)

Maximize(g(x), cons, assume = nonnegative, initialpoint = ip)

[16.4237121534987303, [x[1] = HFloat(0.5028668384263846), x[2] = HFloat(0.49713316157361515), x[3] = HFloat(1.4971331615736156), x[4] = HFloat(3.62694858829851e7)]]

(8)

Maximize(g(x), {x[1]+x[2] <= 1, x[1]+x[3] <= 2}, assume = nonnegative, initialpoint = ip)

[16.4237121534987303, [x[1] = HFloat(0.5028668384263846), x[2] = HFloat(0.49713316157361515), x[3] = HFloat(1.4971331615736156), x[4] = HFloat(3.62694858829851e7)]]

(9)

 

Maximize(sum(ln(x[i]), i = 1 .. 4), {x[1]+x[2] <= 1, x[1]+x[3] <= 2}, assume = nonnegative, initialpoint = ip)

[16.4237121534987303, [x[1] = HFloat(0.5028668384263846), x[2] = HFloat(0.49713316157361515), x[3] = HFloat(1.4971331615736156), x[4] = HFloat(3.62694858829851e7)]]

(10)

``

NULL

Download opt.mw

ranges and scaling?? See the attached.

BTW the gridlines do not appear in an actual Maple worksheet, so the plot *looks* much nicer

restart;
TWeq:=2.96736996560705*10^(-12)*p^2+1.31319840299485*10^(-13)*t^2-8.89549693662593*10^(-7)*p+8.53128393394231*10^(-7)*t-3.65558815509970*10^(-30)*p*t-1 = 0;
plots:-implicitplot(TWeq, p = -8*10^6..2*10^6, t=-8*10^6..2*10^6, scaling=constrained)

0.2967369966e-11*p^2+0.1313198403e-12*t^2-0.8895496937e-6*p+0.8531283934e-6*t-0.3655588155e-29*p*t-1 = 0

 

 

 

 

Download ell.mw

as in the attached

restart; K1 := 10^(-2.12); K2 := 10^(-7.21); K3 := 10^(-12.32); Kw := 10^(-14); NaCl := .5004/(58.44); NaH2PO4 := .1092/(119.98); Na2HPO4 := 1.214/(141.96); Why := Kw = H*OH; Do := K1 = H*H2PO4/H3PO4; It := K2 = H*HPO4/H2PO4; Today := K3 = H*PO4/HPO4; When := Cl = NaCl; Its := Na = NaH2PO4+2*Na2HPO4+NaCl; Due := Na+H = Cl+3*PO4+2*HPO4+H2PO4+OH; Tomorrow := Na2HPO4+NaH2PO4 = H3PO4+H2PO4+HPO4+PO4; z := fsolve([Why, Do, It, Today, When, Its, Due, Tomorrow]); eval(H2PO4, z)

{Cl = 0.8562628337e-2, H = -0.2274751800e-1, H2PO4 = -0.4733982564e-2, H3PO4 = 0.1419582613e-1, HPO4 = 0.1283194935e-7, Na = 0.2657618944e-1, OH = -0.4396084003e-12, PO4 = -0.2699968014e-18}

 

-0.4733982564e-2

(1)

 

Download ChemNum.mw

and I'm not sure what it is!

Luckily, you are only using it to produce the vertical plane x=y, which is trivial to do uing the plottools:polygon() command. Everything now works in the attached. I did change some colors and radii on the tubeplots, to make clearer (to me) what was going on. You may want to change these back, or reset to other values. Same issue with "scaling", Some of the plots are on rather different scales, which makes combining them look a bit odd - you may want to consider some different plot ranges, or even scaling=constrained.

The actual plots look a lot better in the worksheet than they render on this site - honest!
 

Relud surface(2)

 

 

restart;
with(plots):
with(plottools):

G1 := tubeplot([4*t, -1, t], t = -10 .. 10, radius = 2, color = blue):
G2 := tubeplot([1, t, t^2], t = -15 .. 15, radius = 2, color = cyan):
G3 := plottools:-polygon([[-5, -5, -4], [-5, -5, 4], [5, 5, 4], [5, 5, -4]], color = red):
bk := display([G1, G2, G3]);

 

Line:=(q,t)->[4*q+(1-4*q)*t,-1+(1-4*q)*t,q+(16*q^2-q)*t]:
F:= Q-> display(  line(Line(Q,-1.5),Line(Q,+2.5),color=red,thickness=5),
                  plot3d([4*q+y+1,y,q+(y+1)*(16*q^2-q)/(1-4*q)], q=-4..Q, y=(-2.5+6*q)..(1.5-10*q),
                         style=patchcontour, lightmodel=light4, axes=framed)
                        ):

animate(F, [q], q=-4..4, frames=50, lightmodel=light4,background=bk);

 

NULL

Download oddPlot.mw

something like this?

(which won't display "inline", becauase this site is broken AGAIN

fit.mw

The attached shows three of them

  1. using the Student[MultivariateCalculus]() package;
  2. using the LinearAlgebra() package
  3. using the geom3d() package

(And this site is failing to produce an inline display of wroksheets -AGAIN!)

ptline.mw

 

the procedures f1() and f2() in generating equations in e1.. e6, however you are calling these with "square" brackets, eg f1[n]. These should be "round" brackets. In Maple "square" brackets are only used to access entries in an indexable data container

I have fixed this problem in the attached, where I have also restricted the extent of the main "for" loop, just to avoid generating too much output for this page,. Also, just for display purposes, I slightly modified the format strings in the main printf() statement, just for clarity.

For some reason this worksheet won't display inline here, so will you have to download to check

fsol.mw

using PDEtools:-dchange(), as shown in the attached?

restart

"T(x,eta):=Te+Te*(1-g(x,eta))*phi(x);  PDE:=diff(T(x,eta),eta$2)+ f(x,eta)*diff(T(x,eta),eta)=x*(diff(f(x,eta),eta)*diff(T(x,eta),x)-diff(T(x,eta),eta)*diff(f(x,eta),x))    "

proc (x, eta) options operator, arrow, function_assign; Te+Te*(1-g(x, eta))*phi(x) end proc

 

-Te*(diff(diff(g(x, eta), eta), eta))*phi(x)-f(x, eta)*Te*(diff(g(x, eta), eta))*phi(x) = x*((diff(f(x, eta), eta))*(-Te*(diff(g(x, eta), x))*phi(x)+Te*(1-g(x, eta))*(diff(phi(x), x)))+Te*(diff(g(x, eta), eta))*phi(x)*(diff(f(x, eta), x)))

(1)

NULL

NULL

NULL

PDEtools:-dchange(eta = (ue/(c*x))^(1/2)*y, PDE, [y], params = [ue, c])

-Te*(diff(diff(g(y, x), y), y))*c*x*phi(x)/ue-f(y, x)*Te*(diff(g(y, x), y))*phi(x)/(ue/(c*x))^(1/2) = x*((diff(f(y, x), y))*(-Te*(diff(g(y, x), x))*phi(x)+Te*(1-g(y, x))*(diff(phi(x), x)))/(ue/(c*x))^(1/2)+Te*(diff(g(y, x), y))*phi(x)*(diff(f(y, x), x))/(ue/(c*x))^(1/2))

(2)

NULL

Download chngeVar.mw

you could use complexplot(), as in the attached

restart

with(plottools)

with(plots)

with(CurveFitting)

Digits := 10

with(GaussInt)

w := GInearest(0+I)

I

(1)

NULL

"f(t):=7.0*(e)^((-(t-13180)^(2))/(2000000))+4.7*(e)^((-(t-16000)^(2))/(3200000)):"

p1 := plot(f(t), t = 0 .. 20000, color = green); plots[display]({p1})

 

NULL

D1 := 15

epsilon := 200000

L := 6500

v := .7

.7

(2)

t := 10000

10000

(3)

i = sqrt(-1)

i = I

(4)

"k(n) := evalf((2 *Pi*n)/(L))"

proc (n) options operator, arrow, function_assign; evalf(2*Pi*n/L) end proc

(5)

f(n) = (int(f(t)*exp(-w*k(n)*x), x = 0 .. L))/L

``

"C(x, t) :=  (&sum;) exp(-v* t*k(n)- D1 *t*(k(n)^())^(2)- epsilon *t*(k(n))^(4)) *f(n)* exp(w*k(n)* x)"

proc (x, t) options operator, arrow, function_assign; sum(exp(-v*t*k(n)-D1*t*k(n)^2-varepsilon*t*k(n)^4)*f(n)*exp(w*k(n)*x), n = 1 .. 10) end proc

(6)

uu10000 := [seq(evalf(C(L-j, t)), j = 0 .. 6500, 100)]; complexplot(uu10000, x = min(`~`[Re](uu10000)) .. max(`~`[Re](uu10000)), scaling = constrained, style = point)

 

``

Download cplxplt.mw

is to solve the ODE given as f'''+ff''-f'^2-Mf'=0, with f(0)=0, f'(0)=1, and f'(5)=0 with M=0.5

So let's just solve it, as shown in the attached.

  restart:
  with(plots):

  M:=0.5;
  ode:= diff(f(x), x$3)+f(x)*diff(f(x), x$2)-diff(f(x),x)^2-M*diff(f(x),x)=0;
  conds:= f(0)=0, D(f)(0)=1, D(f)(5 )-0;
  ans:=dsolve( [ode, conds], numeric);
  odeplot(ans, [x, f(x)], x=0..5);
  

.5

 

diff(diff(diff(f(x), x), x), x)+f(x)*(diff(diff(f(x), x), x))-(diff(f(x), x))^2-.5*(diff(f(x), x)) = 0

 

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

 

proc (x_bvp) local res, data, solnproc, _ndsol, outpoint, i; option `Copyright (c) 2000 by Waterloo Maple Inc. All rights reserved.`; _EnvDSNumericSaveDigits := Digits; Digits := 15; if _EnvInFsolve = true then outpoint := evalf[_EnvDSNumericSaveDigits](x_bvp) else outpoint := evalf(x_bvp) end if; data := Array(1..4, {(1) = proc (outpoint) local X, Y, YP, yout, errproc, L, V, i; option `Copyright (c) 2000 by Waterloo Maple Inc. All rights reserved.`; X := Vector(13, {(1) = .0, (2) = .3648083980208248, (3) = .739453123146836, (4) = 1.1294231031240023, (5) = 1.5404396951722952, (6) = 1.9674194416734547, (7) = 2.4067451277606597, (8) = 2.850771504517652, (9) = 3.2971678108063687, (10) = 3.7442923732839364, (11) = 4.191663594165653, (12) = 4.620140941040035, (13) = 5.0}, datatype = float[8], order = C_order); Y := Matrix(13, 3, {(1, 1) = .0, (1, 2) = 1.0, (1, 3) = -1.2249064035636923, (2, 1) = .294194453432655, (2, 2) = .6396128765313516, (2, 3) = -.7836115514965457, (3, 1) = .48635488180500724, (3, 2) = .40415085049570865, (3, 3) = -.49534834686074514, (4, 1) = .6116367675260377, (4, 2) = .250541601628673, (4, 3) = -.3073631245893844, (5, 1) = .6925060924519302, (5, 2) = .15125045647457197, (5, 3) = -.18593456512252354, (6, 1) = .7427378377645504, (6, 2) = 0.893905555320538e-1, (6, 3) = -.110380533553483, (7, 1) = .7730294979143654, (7, 2) = 0.51843487043024994e-1, (7, 3) = -0.6463752743281678e-1, (8, 1) = .7906832855559992, (8, 2) = 0.29654107059807183e-1, (8, 3) = -0.3774106284132804e-1, (9, 1) = .8007523087394047, (9, 2) = 0.1661644456196504e-1, (9, 3) = -0.22098378262475696e-1, (10, 1) = .8063158313666053, (10, 2) = 0.8938405451087109e-2, (10, 3) = -0.1307531449846363e-1, (11, 1) = .8092030889605896, (11, 2) = 0.4352616339821983e-2, (11, 3) = -0.7907830531921609e-2, (12, 1) = .810440946160643, (12, 2) = 0.1626507257535621e-2, (12, 3) = -0.5077258781403503e-2, (13, 1) = .8107323262402987, (13, 2) = .0, (13, 3) = -0.3612766608982084e-2}, datatype = float[8], order = C_order); YP := Matrix(13, 3, {(1, 1) = 1.0, (1, 2) = -1.2249064035636923, (1, 3) = 1.5, (2, 1) = .6396128765313516, (2, 2) = -.7836115514965457, (2, 3) = .9594452421864268, (3, 1) = .40415085049570865, (3, 2) = -.49534834686074514, (3, 3) = .6063284218940224, (4, 1) = .250541601628673, (4, 2) = -.3073631245893844, (4, 3) = .37603648294155106, (5, 1) = .15125045647457197, (5, 2) = -.18593456512252354, (5, 3) = .2272627479658001, (6, 1) = 0.893905555320538e-1, (6, 2) = -.110380533553483, (6, 3) = .13466974800716747, (7, 1) = 0.51843487043024994e-1, (7, 2) = -0.6463752743281678e-1, (7, 3) = 0.7857620604810917e-1, (8, 1) = 0.29654107059807183e-1, (8, 2) = -0.3774106284132804e-1, (8, 3) = 0.4554764716317479e-1, (9, 1) = 0.1661644456196504e-1, (9, 2) = -0.22098378262475696e-1, (9, 3) = 0.26279655923937467e-1, (10, 1) = 0.8938405451087109e-2, (10, 2) = -0.1307531449846363e-1, (10, 3) = 0.1509193089776011e-1, (11, 1) = 0.4352616339821983e-2, (11, 2) = -0.7907830531921609e-2, (11, 3) = 0.8594294332320503e-2, (12, 1) = 0.1626507257535621e-2, (12, 2) = -0.5077258781403503e-2, (12, 3) = 0.4930717565329715e-2, (13, 1) = .0, (13, 2) = -0.3612766608982084e-2, (13, 3) = 0.29289866770633205e-2}, datatype = float[8], order = C_order); errproc := proc (x_bvp) local outpoint, X, Y, yout, L, V, i; option `Copyright (c) 2000 by Waterloo Maple Inc. All rights reserved.`; Digits := 15; outpoint := evalf(x_bvp); X := Vector(13, {(1) = .0, (2) = .3648083980208248, (3) = .739453123146836, (4) = 1.1294231031240023, (5) = 1.5404396951722952, (6) = 1.9674194416734547, (7) = 2.4067451277606597, (8) = 2.850771504517652, (9) = 3.2971678108063687, (10) = 3.7442923732839364, (11) = 4.191663594165653, (12) = 4.620140941040035, (13) = 5.0}, datatype = float[8], order = C_order); Y := Matrix(13, 3, {(1, 1) = .0, (1, 2) = .0, (1, 3) = 0.1273975608644069e-9, (2, 1) = -0.9699958654485496e-7, (2, 2) = 0.11898009711795105e-6, (2, 3) = -0.14531470840848373e-6, (3, 1) = -0.10320306463506545e-6, (3, 2) = 0.12691811465663893e-6, (3, 3) = -0.154678061248875e-6, (4, 1) = -0.7871067918956964e-7, (4, 2) = 0.9744291713856059e-7, (4, 3) = -0.11818172410834043e-6, (5, 1) = -0.4865437588049475e-7, (5, 2) = 0.613589205139654e-7, (5, 3) = -0.7358044074399424e-7, (6, 1) = -0.21030326071357806e-7, (6, 2) = 0.28453569331208245e-7, (6, 3) = -0.3288511507213022e-7, (7, 1) = -0.12272437930054753e-8, (7, 2) = 0.5323548042522246e-8, (7, 3) = -0.4173055146179112e-8, (8, 1) = 0.9931557233651253e-8, (8, 2) = -0.7035041293820884e-8, (8, 3) = 0.11351363226824974e-7, (9, 1) = 0.14250962440952772e-7, (9, 2) = -0.10817614549677532e-7, (9, 3) = 0.1640579916342742e-7, (10, 1) = 0.14719058046583237e-7, (10, 2) = -0.965367179587185e-8, (10, 3) = 0.15469412664652983e-7, (11, 1) = 0.13753746410766736e-7, (11, 2) = -0.6455085789179595e-8, (11, 3) = 0.12146772829413411e-7, (12, 1) = 0.1268333027160663e-7, (12, 2) = -0.2912324487441674e-8, (12, 3) = 0.8508819311976027e-8, (13, 1) = 0.12136283725856026e-7, (13, 2) = .0, (13, 3) = 0.5688098437791177e-8}, datatype = float[8], order = C_order); if not type(outpoint, 'numeric') then if outpoint = "start" or outpoint = "left" then return X[1] elif outpoint = "right" then return X[13] elif outpoint = "order" then return 6 elif outpoint = "error" then return HFloat(1.54678061248875e-7) elif outpoint = "errorproc" then error "this is already the error procedure" elif outpoint = "rawdata" then return [3, 13, [f(x), diff(f(x), x), diff(diff(f(x), x), x)], X, Y] else return ('procname')(x_bvp) end if end if; if outpoint < X[1] or X[13] < outpoint then error "solution is only defined in the range %1..%2", X[1], X[13] end if; V := array([1 = 4, 2 = 0]); if Digits <= trunc(evalhf(Digits)) then L := Vector(4, 'datatype' = 'float'[8]); yout := Vector(3, 'datatype' = 'float'[8]); evalhf(`dsolve/numeric/lagrange`(13, 3, X, Y, outpoint, var(yout), var(L), var(V))) else L := Vector(4, 'datatype' = 'sfloat'); yout := Vector(3, 'datatype' = 'sfloat'); `dsolve/numeric/lagrange`(13, 3, X, Y, outpoint, yout, L, V) end if; [x = outpoint, seq('[f(x), diff(f(x), x), diff(diff(f(x), x), x)]'[i] = yout[i], i = 1 .. 3)] end proc; if not type(outpoint, 'numeric') then if outpoint = "start" or outpoint = "left" then return X[1] elif outpoint = "method" then return "bvp" elif outpoint = "right" then return X[13] elif outpoint = "order" then return 6 elif outpoint = "error" then return HFloat(1.54678061248875e-7) elif outpoint = "errorproc" then return eval(errproc) elif outpoint = "rawdata" then return [3, 13, "depnames", X, Y, YP] else error "non-numeric value" end if end if; if outpoint < X[1] or X[13] < outpoint then error "solution is only defined in the range %1..%2", X[1], X[13] end if; if Digits <= trunc(evalhf(Digits)) and (_EnvInFsolve <> true or _EnvDSNumericSaveDigits <= trunc(evalhf(Digits))) then V := array( 1 .. 6, [( 1 ) = (7), ( 2 ) = (0), ( 3 ) = (false), ( 4 ) = (false), ( 5 ) = (false), ( 6 ) = (false)  ] ); L := Matrix(7, 2, {(1, 1) = .0, (1, 2) = .0, (2, 1) = .0, (2, 2) = .0, (3, 1) = .0, (3, 2) = .0, (4, 1) = .0, (4, 2) = .0, (5, 1) = .0, (5, 2) = .0, (6, 1) = .0, (6, 2) = .0, (7, 1) = .0, (7, 2) = .0}, datatype = float[8], order = C_order); yout := Vector(3, {(1) = .0, (2) = .0, (3) = .0}, datatype = float[8]); evalhf(`dsolve/numeric/hermite`(13, 3, X, Y, YP, outpoint, var(yout), var(L), var(V))) else if _EnvInFsolve = true then Digits := _EnvDSNumericSaveDigits end if; V := array( 1 .. 6, [( 1 ) = (7), ( 2 ) = (0), ( 3 ) = (false), ( 4 ) = (false), ( 5 ) = (false), ( 6 ) = (false)  ] ); L := Matrix(7, 2, {(1, 1) = 0., (1, 2) = 0., (2, 1) = 0., (2, 2) = 0., (3, 1) = 0., (3, 2) = 0., (4, 1) = 0., (4, 2) = 0., (5, 1) = 0., (5, 2) = 0., (6, 1) = 0., (6, 2) = 0., (7, 1) = 0., (7, 2) = 0.}, order = C_order); yout := Vector(3, {(1) = 0., (2) = 0., (3) = 0.}); `dsolve/numeric/hermite`(13, 3, X, Y, YP, outpoint, yout, L, V) end if; [outpoint, seq(yout[i], i = 1 .. 3)] end proc, (2) = Array(0..0, {}), (3) = [x, f(x), diff(f(x), x), diff(diff(f(x), x), x)], (4) = 0}); solnproc := data[1]; if not type(outpoint, 'numeric') then if outpoint = "solnprocedure" then return eval(solnproc) elif member(outpoint, ["start", "left", "right", "errorproc", "rawdata", "order", "error"]) then return solnproc(x_bvp) elif outpoint = "sysvars" then return data[3] elif procname <> unknown then return ('procname')(x_bvp) else _ndsol := pointto(data[2][0]); return ('_ndsol')(x_bvp) end if end if; try res := solnproc(outpoint); [x = res[1], seq('[f(x), diff(f(x), x), diff(diff(f(x), x), x)]'[i] = res[i+1], i = 1 .. 3)] catch: error  end try end proc

 

 

 

Download solODE.mw

which if you differentiate it, will give a 7-th order polynomial. Setting the latter equal to zero will rsult in seven roots (ie seven stationary points), although as mmcdara has pointed out three of these occur at x=-1. \Whether you consider this a single stationary point (with a multiplicity of 3), or three coincident stationary points is a quibble about language

The other four roots are complex

Check the tickmarks!

The problem appears to be that you have two sets of horizontal gridlines - one for each of the vertical axes in the dualaxis plot. Since the left and right tickmarks do not align, these two sets of gridlines are interleaved, giving an "illusion" of a weird vertical scale.

In a dualaxis plot, I can't think of any way to make the left and right tickmarks "align", since they are likely to be on very different scales.

Even if one could come up with a method of having horizontal gridlines only for the left axis (and I don't know if this can be done!), it means that these gridlines would not coincide with the tickmarks on the right axis - which would make the plot look a bit weird

wrong with the fsolve() command using the 'complex' option?

https://www.mapleprimes.com/questions/235892-How-Do-I-Solve-A-Differential-Equation-In-Maple

3 4 5 6 7 8 9 Last Page 5 of 207