I do not think the polynomial you specify, with all coefficients positive, has a positive root.  Assuming positive does not include zero or complex.

restart
assume(a > 0, b > 0, c > 0, b^2 >= 4*a*c)
eqn := a*x^2+b*x+c = 0
soln := solve(eqn, x)
for s in soln do
  if verify(s, 0, greater_than) then
    print("Greater than", s)
  else
    print("Negative", s)
  end if
end do

@Preben Alsholm The NASA paper,

https://ntrs.nasa.gov/archive/nasa/casi.ntrs.nasa.gov/19690021375.pdf

was aimed at deriving constants for use in RK methods.  How the 'exact solution' was determined was not mentioned.  The problem I was looking at is given on page 38. 

A second PDE example if given on page 40 of the paper.  But, pdetest does not verify that this is a solution...

Tom Dean

Use 'create a new document block or enclose the current section in one' followed by 'insert a maple prompt after the current execution group'

Sorry for the noise.

Tom Dean

@Carl Love I am on x64 Linux. (core i7, OC 4.2GHz, 16G RAM)

I have default axes names defined in .mapleinit.

I use emacs as my interface to command line maple.  So, the double underscore remains a double underscore.  With xmaple,  the double underscore produces subscript, visually.  This should be a keyboard shortcut.  But, in xmaple, using

p__1 - p[1];

does not result in zero.  I wonder why, since they are visually the same.

p__1 := 7; p[1];

Shows that things that are visually the same are not.

This is the same with command line maple.

@Carl Love Old UNIX habit...

(a,b) := op(calcView(relBrgLinePlot,relTgtPositPlot,relTgtLinePlot));
                    a, b := 0. .. 95., -61.58179461 .. 16.

viewArg := [floor(op(1,a))-10..ceil(op(2,a))+10,
            floor(op(1,b))-10..ceil(op(2,b))+10];
                      viewArg := [-10 .. 105, -72 .. 26]

display([relBrgLinePlot,relTgtPositPlot,relTgtLinePlot],
        title="Relative Motion",
        labels=["Miles, +North","Miles +East"],
        labeldirections=["horizontal","vertical"],
        gridlines = true,
        axes=boxed,
        view=viewArg);

calcView := proc(itmList)
local v, x, y, itm, minx, maxx, miny, maxy;
local absminx, absmaxx, absminy, absmaxy, thisitm;
    absminx := 10000;
    absmaxx := -10000;
    absminy := 10000;
    absmaxy := -10000;
    for itm in itmList do
        v := op(-1,itm):
        x := op(1,v):
        y := op(2,v):
        (minx, maxx, miny, maxy) := (op(1,x), op(2,x), op(1,y), op(2,y));
        ##print(minx, maxx, miny, maxy):
        absminx := min(minx, absminx);
        absmaxx := max(maxx, absmaxx);
        absminy := min(miny, absminy);
        absmaxy := max(maxy, absmaxy);
    end do;
    [floor(absminx)-10, ceil(absmaxx)+10,\
     floor(absminy)-10, ceil(absmaxy)+10]:
end proc:

> calcView([relBrgLinePlot,relTgtPositPlot,relTgtLinePlot]);
                              [-10, 105, -72, 26]

> calcView([geoSensPositPlot,geoSensLinePlot,\
            geoBrgLinePlot,geoTgtPositPlot,geoTgtLinePLot]);
                              [-10, 105, -10, 98]

@Mac Dude Thanks, I use Maple2016 on Linux x64, with plotsetup(x11), no gridlines.

Tom Dean

@Markiyan Hirnyk With maple2016, I get equal values, at least to the number of digits displayed.

J:= Int(I*sqrt((R*exp(I*theta)+1)/(R*exp(I*theta)-a)), theta= 0..Pi) ;

> K := evalf(eval(J, [R = 1/2, a = 1/4]));                                     
                       K := 2.677293769 + 3.506183701 I

> int(eval(I*sqrt((R*exp(I*theta)+1)/(R*exp(I*theta)-a)), [R = .5, a = .25]), theta = 0 .. Pi,numeric);
                          2.677293769 + 3.506183701 I

Tom Dean

@Carl Love Thanks. Your first post on the plot m vs y thread was better yet because of more hardware use.

@Carl Love Thanks, I learned some things.

@Carl Love I see two ways of the package function long form.  I think most help uses [].

I looked at Statistics and got interested, but, maybe confused?

restart;
with(Statistics):
y := int(1/(-0.4016e-1*m^(2/3)-0.211e-3*m^(5/3)), m);
## fix complex values
f := (x) -> abs(evalf(subs(m=x,y)));

n := 500;

## mean 1, std dev 0.5
X := RandomVariable(Normal(1, .5)):
M := Sample(X,n):
Y := [seq(f(M[idx]),idx=1..n)]:

A:=[seq(idx,idx=1..n)]:
plot(A,M,style=point,symbol= circle, symbolsize= 1,title="Values of M");
plot(A,Y,style=point,symbol= circle, symbolsize= 1,title="Values of Y");
plot(Y,M,style=point,symbol= circle, symbolsize= 1,title="M vs Y");
DataSummary(M);DataSummary(Y);
Data:=Matrix(n,2,[seq([Y[idx],M[idx]][],idx=1..n)]);

plot(
     [[y, m, m= 0..2], Data], labels= ['y', m],
     style= [line, point], color= [red, black],
     symbol= point, symbolsize= 1
);

@Carl Love What did you use to generate the data.  I have been looking at the Statistics package for an hour, with little success.

Tom Dean

@Carl Love Thanks, I understand that.  My difficulty was the sensitivity of the numerical result
to relatively minor changes to the input values.

## bearing.mpl, solve the target motion problem with bearings only.
##
## Consider a sensor platform moving through points (x,y) at times
## t[1..4] with the target bearings, Brg[1..4] taken at times t[1..4]
## with the target proceeding along a constant course and speed.
##
## time t, bearing line slope m, sensor position (x,y) are known
## values.
##
## Since this is a generated problem the target position at time t is
## provided to compare with the results.
##
#########################################################################
##
restart;
##
genKnownValues := proc()
    description "set the known values",
    "t - relative time",
    "x - sensor x location at time t[i]",
    "y - sensor y location at time t[i]",
    "m - slope of the bearing lines at time t[i]",
    "tgtPosit - target position at time t[i]";
    global t, mBrg, mCalc, x, y, tgtPosit;
    local dt, Cse, Spd, Brg, A, B, C, R, X;
    local tgtX, tgtY, tgtRange, tgtCse, tgtSpd;
## relative and delta time
    t := [0, 1+1/2, 3, 3+1/2];
    dt := [0, seq(t[idx]-t[idx-1],idx=2..4)];
## sensor motion
    Cse := [90, 90, 90, 50] *~ Pi/180; ## true heading
    Spd := [15, 15, 15, 22];  ## knots
## bearings to the target at time t
    Brg := [10, 358, 340, 330] *~ (Pi/180);
## calculate the sensor position vs time
    x := ListTools[PartialSums](dt *~ Spd *~ map(cos, Cse));
    y := ListTools[PartialSums](dt *~ Spd *~ map(sin, Cse));
## target values  start the target at a known (x,y) position at a
## constant course and speed
    tgtRange := 95+25/32; ## miles at t1, match octave value...
    tgtCse := 170 * Pi/180; ## course
    tgtSpd := 10; ## knots
    ## tgtX := tgtRange*cos(Brg[1]);
    tgtX := 94 + 1/(3 + 1/(3 + 1/(2 + 1/3)));
    tgtX := tgtX +~ ListTools[PartialSums](dt *~ tgtSpd *~ cos(tgtCse));
    ## tgtY := tgtRange*sin(Brg[1]);
    tgtY := 17 + 1/(-4 + 1/(-3 + 1/10));
    tgtY := tgtY +~ ListTools[PartialSums](dt *~ tgtSpd *~ sin(tgtCse));
## return target position vs time as a matrix
    tgtPosit:=Matrix(4,2,[seq([tgtX[idx],tgtY[idx]][],idx=1..4)]);
## slope of the bearing lines
    mBrg := map(tan,Brg);
    mCalc := [seq((tgtY[idx]-y[idx])/(tgtX[idx] - x[idx]),idx=1..4)];
end proc:
##
#########################################################################
## t[], m[], x[], and y[] are known values
##
## equation of the bearing lines
eq1 := tgtY[1] - y[1]    = m[1]*(tgtX[1]-x[1]):
eq2 := tgtY[2] - y[2]    = m[2]*(tgtX[2]-x[2]):
eq3 := tgtY[3] - y[3]    = m[3]*(tgtX[3]-x[3]):
eq4 := tgtY[4] - y[4]    = m[4]*(tgtX[4]-x[4]):
## target X motion along the target line
eq5 := tgtX[2] - tgtX[1] = tgtVx*(t[2]-t[1]):
eq6 := tgtX[3] - tgtX[2] = tgtVx*(t[3]-t[2]):
eq7 := tgtX[4] - tgtX[3] = tgtVx*(t[4]-t[3]):
## target Y motion along the target line
eq8 := tgtY[2] - tgtY[1] = tgtVy*(t[2]-t[1]):
eq9 := tgtY[3] - tgtY[2] = tgtVy*(t[3]-t[2]):
eq10:= tgtY[4] - tgtY[3] = tgtVy*(t[4]-t[3]):
##
#########################################################################
##
## solve the equations
eqs  := {eq1,eq2,eq3,eq4,eq5,eq6,eq7,eq8,eq9,eq10}:

Sol:= solve(eqs, {tgtVx, tgtVy, seq([tgtX[k], tgtY[k]][], k= 1..4)}):
##
SolSet := [tgtX[4]=eval(tgtX[4], Sol),
           tgtY[4]=eval(tgtY[4], Sol),
           tgtVx=eval(tgtVx, Sol),
           tgtVy=eval(tgtVy, Sol)]:

genKnownValues():
with(geometry):
m:= mBrg:
point(P1,[eval(tgtX[4],SolSet),eval(tgtY[4],SolSet)]):
m := mCalc:
point(P2,[eval(tgtX[4],SolSet),eval(tgtY[4],SolSet)]):
evalf(distance(P1,P2));

Is there a way to select the product without reposting the question?

I needed to reboot after installing the xfonts...