Advanced Search

dharr

Dr. David Harrington

9047 Reputation

22 Badges

21 years, 200 days
University of Victoria
Professor or university staff
Victoria, British Columbia, Canada

Social Networks and Content at Maplesoft.com

Maple Application Center
Contact Dr. David Harrington
I am a retired professor of chemistry at the University of Victoria, BC, Canada. My research areas are electrochemistry and surface science. I have been a user of Maple since about 1990.

MaplePrimes Activity


These are replies submitted by dharr

File formats...

dharr 9047
5 hours ago
1 0

These are not files listed on the File format help page, so you are out of luck for an easy way. I guess if you know the details of the binary file you could write a routine to read them, but I don't think that would count as "easy". 

didn't see this...

dharr 9047
Yesterday at 5:46 PM
0 0

Just installed 2026.0 in windows 11; run as administrator to install; didn't import previous preferences. All is working fine.

comment...

dharr 9047
Yesterday at 1:26 PM
0 0

@KIRAN SAJJAN A quick look suggests it is reasonable, assuming the R1(z) and R2(z) functions are correct. I believe you said the BCs depended on z, but I do not see that dependence, though it might be implicit.

confusion...

dharr 9047
March 11 2026
0 0

@Suryakanth @KIRAN SAJJAN You set up the odes for other conditions, so why not set up for different R or z and solve them? Sorry, I don't understand what you want.

z dependence...

dharr 9047
March 11 2026
0 0

@KIRAN SAJJAN Sorry, the equations you provided do not contain z, so I'm confused. Do you mean it is another independent variable, and you want to solve a PDE? (I don't want to read and figure out the paper.)

fixed in 2025.2...

dharr 9047
March 11 2026
0 0

@WD0HHU This seemed to be fixed in 2025.2

confusing with parametric plot...

dharr 9047
March 11 2026
0 0

@salim-barzani Use a set for the 3 plots otherwise it will think it is a parametric plot.

Dr.D-N.mw

one way...

dharr 9047
March 11 2026
0 0

You can directly use the procedure for w(eta) returned by dsolve to evaluate w(0.999) etc and the integrals.

parentartery_and_daughter_artery_error.mw

numerical test...

dharr 9047
March 11 2026
0 0

@salim-barzani That explains why odetest didn't work; the rest is OK. The sum of kinetic energy and potential energy should be approximately constant, not each separately.

You undoubtedly know some better numerical values for the parameters. Here the variation in energy is about 20%. You could also just plot the approximate solution and the exact solution (which you know?) together and see how they look.

> 

restart
> 

with(PDEtools)
> 

undeclare(prime, quiet)
> 

declare(u(x, t), quiet); declare(U(xi), quiet); declare(V(xi), quiet); declare(Omega(x, y, t), quiet); declare(Q(x, t), quiet); declare(R(x, t), quiet); declare(CHI(x, t), quiet); declare(Z(x, t), quiet)
> 

pde1 := I*(diff(Omega(x, y, t), t))-I*(diff(Omega(x, y, t), x))+alpha[1]*(diff(Omega(x, y, t), `$`(x, 2)))+alpha[2]*(diff(Omega(x, y, t), `$`(y, 2)))+alpha[3]*(diff(Omega(x, y, t), x, y))-alpha[4]*Omega(x, y, t)*U(p*x+q*y-t*v)^2
> 

G1 := U(p*x+q*y-t*v) = U(xi); G2 := (D(U))(p*x+q*y-t*v) = diff(U(xi), xi); G3 := ((D@@2)(U))(p*x+q*y-t*v) = diff(U(xi), `$`(xi, 2)); G4 := ((D@@3)(U))(p*x+q*y-t*v) = diff(U(xi), `$`(xi, 3)); G5 := ((D@@4)(U))(p*x+q*y-t*v) = diff(U(xi), `$`(xi, 4)); G6 := ((D@@5)(U))(p*x+q*y-t*v) = diff(U(xi), `$`(xi, 5))
> 

T := xi = p*x+q*y-t*v; T1 := Omega(x, y, t) = U(p*x+q*y-t*v)*exp(I*(-r*t+theta*y+w*x+phi))
> 

P1 := I*(diff(Omega(x, y, t), t))-I*(diff(Omega(x, y, t), x))+alpha[1]*(diff(Omega(x, y, t), `$`(x, 2)))+alpha[2]*(diff(Omega(x, y, t), `$`(y, 2)))+alpha[3]*(diff(Omega(x, y, t), x, y))-alpha[4]*Omega(x, y, t)*U(p*x+q*y-t*v)^2
> 

P11 := eval(P1, {T, T1})
> 

P111 := subs({G1, G2, G3, G4, G5, G6}, P11)
> 

pde1 := P111 = 0
> 

D1 := numer(lhs(pde1))*denom(rhs(pde1)) = numer(rhs(pde1))*denom(lhs(pde1)); simplify(%, 'symbolic'); ode := collect(%, {U(xi), diff(diff(U(xi), xi), xi)})

Eq. 37

> 

ode1 := expand(lhs(ode)/exp(I*(-r*t+theta*y+w*x+phi)))
> 

ode11 := evalc(Re(ode1))
> 

ode111 := collect(ode11, {U(xi), diff(diff(U(xi), xi), xi)})

-U(xi)^3*alpha[4]+(-theta^2*alpha[2]-theta*w*alpha[3]-w^2*alpha[1]+r+w)*U(xi)+(p^2*alpha[1]+p*q*alpha[3]+q^2*alpha[2])*(diff(diff(U(xi), xi), xi))

Eq, 5.1

> 

ode1111 := collect(ode111/DEtools:-dcoeffs(ode111, U(xi))[1], U(xi))

Note that each term in the Hamiltonian integrand when differentiated wrt U(xi) gives a term in ode. (I have no idea why we do this.) So differentiating wrt xi would give a term in the ode multiplied by diff(U(xi),xi) (using the chain rule). So we can set up a differential equation for the integrand f(xi)

> 

def := diff(f(xi), xi) = ode1111*(diff(U(xi), xi))

We want the first integral of this

> 

DEtools:-firint(def, f(xi)); integrand := collect(solve(eval(%, {_C1 = 0, c__1 = 0}), f(xi)), diff, simplify)

From which we can extract K and P.

> 

K := op(1, integrand); P := collect(-op(2, integrand), U(xi))

Eq. 5.6

> 

H := K+P

(1/2)*(diff(U(xi), xi))^2+alpha[4]*U(xi)^4/(4*p^2*alpha[1]+4*p*q*alpha[3]+4*q^2*alpha[2])+(2*theta^2*alpha[2]+2*theta*w*alpha[3]+2*w^2*alpha[1]-2*r-2*w)*U(xi)^2/(4*p^2*alpha[1]+4*p*q*alpha[3]+4*q^2*alpha[2])

Eq. 5.7

> 

sol1 := U(xi) = m*cos(delta*xi)

U(xi) = m*cos(delta*xi)

Eq. 5.10

> 

eval(convert(H, D), xi = 0); H0 := collect(eval(%, {U(0) = m, (D(U))(0) = 0}), m)

Eq. 5.11 - looks correct

> 

H1 := eval(H, sol1)-H0

Eq. 5.12 - looks correct

> 

H2 := collect(algsubs(delta*xi = (1/4)*Pi, H1), delta, proc (x) options operator, arrow; collect(x, m) end proc)

Eq. 5.13

> 

delta = sort([solve(H2, delta)])[2]; eqphi := `assuming`([simplify(%)], [p^2*alpha[1]+p*q*alpha[3]+q^2*alpha[2] > 0])

delta = (1/2)*(4*alpha[1]*w^2+(4*theta*alpha[3]-4)*w+3*alpha[4]*m^2+4*theta^2*alpha[2]-4*r)^(1/2)/(p^2*alpha[1]+p*q*alpha[3]+q^2*alpha[2])^(1/2)

> 

sol2 := simplify(eval(sol1, eqphi))

U(xi) = m*cos((1/2)*(4*alpha[1]*w^2+(4*theta*alpha[3]-4)*w+3*alpha[4]*m^2+4*theta^2*alpha[2]-4*r)^(1/2)*xi/(p^2*alpha[1]+p*q*alpha[3]+q^2*alpha[2])^(1/2))

So the kinetic, potential and total energies are

> 

Kapprox := simplify(eval(K, sol2)); Papprox := simplify(eval(P, sol2)); Happrox := simplify(eval(H, sol2))

Which we expect to be approximately constant. Test with some sensible numerical values - her I just choose them all 1

> 

params := `~`[`=`](`minus`(indets(Happrox, name), {xi}), 1)

{m = 1, p = 1, q = 1, r = 1, theta = 1, w = 1, alpha[1] = 1, alpha[2] = 1, alpha[3] = 1, alpha[4] = 1}

> 

Knum := eval(Kapprox, params); Pnum := eval(Papprox, params); Hnum := eval(Happrox, params)

7/48-(7/48)*cos((1/3)*3^(1/2)*7^(1/2)*xi)

(1/12)*cos((1/6)*3^(1/2)*7^(1/2)*xi)^2*(cos((1/6)*3^(1/2)*7^(1/2)*xi)^2+2)

(1/12)*cos((1/6)*3^(1/2)*7^(1/2)*xi)^4-(1/8)*cos((1/6)*3^(1/2)*7^(1/2)*xi)^2+7/24

Total energy varies by about 0.05 compared to K and P p-p amplitudes of 0.25, or 20%.

> 

plot([Knum, Pnum, Hnum], xi = 0 .. 5, view = [default, default], color = [red, blue, black])

> 

NULL

Download u-1.mw

Better...

dharr 9047
March 10 2026
0 0

@salim-barzani Most equations look correct but odetest fails - maybe check some assumptions or try some numerical checks, or maybe the theory doesn't work.

> 

restart
> 

with(PDEtools)
> 

undeclare(prime, quiet)
> 

declare(u(x, t), quiet); declare(U(xi), quiet); declare(V(xi), quiet); declare(Omega(x, y, t), quiet); declare(Q(x, t), quiet); declare(R(x, t), quiet); declare(CHI(x, t), quiet); declare(Z(x, t), quiet)
> 

pde1 := I*(diff(Omega(x, y, t), t))-I*(diff(Omega(x, y, t), x))+alpha[1]*(diff(Omega(x, y, t), `$`(x, 2)))+alpha[2]*(diff(Omega(x, y, t), `$`(y, 2)))+alpha[3]*(diff(Omega(x, y, t), x, y))-alpha[4]*Omega(x, y, t)*U(p*x+q*y-t*v)^2

I*(diff(Omega(x, y, t), t))-I*(diff(Omega(x, y, t), x))+alpha[1]*(diff(diff(Omega(x, y, t), x), x))+alpha[2]*(diff(diff(Omega(x, y, t), y), y))+alpha[3]*(diff(diff(Omega(x, y, t), x), y))-alpha[4]*Omega(x, y, t)*U(p*x+q*y-t*v)^2

> 

G1 := U(p*x+q*y-t*v) = U(xi); G2 := (D(U))(p*x+q*y-t*v) = diff(U(xi), xi); G3 := ((D@@2)(U))(p*x+q*y-t*v) = diff(U(xi), `$`(xi, 2)); G4 := ((D@@3)(U))(p*x+q*y-t*v) = diff(U(xi), `$`(xi, 3)); G5 := ((D@@4)(U))(p*x+q*y-t*v) = diff(U(xi), `$`(xi, 4)); G6 := ((D@@5)(U))(p*x+q*y-t*v) = diff(U(xi), `$`(xi, 5))
> 

T := xi = p*x+q*y-t*v; T1 := Omega(x, y, t) = U(p*x+q*y-t*v)*exp(I*(-r*t+theta*y+w*x+phi))
> 

P1 := I*(diff(Omega(x, y, t), t))-I*(diff(Omega(x, y, t), x))+alpha[1]*(diff(Omega(x, y, t), `$`(x, 2)))+alpha[2]*(diff(Omega(x, y, t), `$`(y, 2)))+alpha[3]*(diff(Omega(x, y, t), x, y))-alpha[4]*Omega(x, y, t)*U(p*x+q*y-t*v)^2
> 

P11 := eval(P1, {T, T1})
> 

P111 := subs({G1, G2, G3, G4, G5, G6}, P11)
> 

pde1 := P111 = 0
> 

D1 := numer(lhs(pde1))*denom(rhs(pde1)) = numer(rhs(pde1))*denom(lhs(pde1)); simplify(%, 'symbolic'); ode := collect(%, {U(xi), diff(diff(U(xi), xi), xi)})

-alpha[4]*U(xi)^3*exp(I*(-r*t+theta*y+w*x+phi))+exp(I*(-r*t+theta*y+w*x+phi))*(-alpha[1]*w^2-(theta*alpha[3]-1)*w-theta^2*alpha[2]+r)*U(xi)+exp(I*(-r*t+theta*y+w*x+phi))*(p^2*alpha[1]+p*q*alpha[3]+q^2*alpha[2])*(diff(diff(U(xi), xi), xi))+I*exp(I*(-r*t+theta*y+w*x+phi))*((theta*alpha[3]+2*w*alpha[1]-1)*p+2*q*theta*alpha[2]+q*w*alpha[3]-v)*(diff(U(xi), xi)) = 0

Eq. 37

> 

ode1 := expand(lhs(ode)/exp(I*(-r*t+theta*y+w*x+phi)))

I*(diff(U(xi), xi))*p*theta*alpha[3]+(2*I)*(diff(U(xi), xi))*p*w*alpha[1]+(2*I)*(diff(U(xi), xi))*q*theta*alpha[2]+I*(diff(U(xi), xi))*q*w*alpha[3]-U(xi)^3*alpha[4]-U(xi)*theta^2*alpha[2]-U(xi)*theta*w*alpha[3]-U(xi)*w^2*alpha[1]+(diff(diff(U(xi), xi), xi))*p^2*alpha[1]+(diff(diff(U(xi), xi), xi))*p*q*alpha[3]+(diff(diff(U(xi), xi), xi))*q^2*alpha[2]-I*(diff(U(xi), xi))*p-I*(diff(U(xi), xi))*v+U(xi)*r+U(xi)*w

> 

ode11 := evalc(Re(ode1))*o

-U(xi)^3*alpha[4]-U(xi)*theta^2*alpha[2]-U(xi)*theta*w*alpha[3]-U(xi)*w^2*alpha[1]+(diff(diff(U(xi), xi), xi))*p^2*alpha[1]+(diff(diff(U(xi), xi), xi))*p*q*alpha[3]+(diff(diff(U(xi), xi), xi))*q^2*alpha[2]+U(xi)*r+U(xi)*w

> 

ode111 := collect(ode11, {U(xi), diff(diff(U(xi), xi), xi)})

-U(xi)^3*alpha[4]+(-theta^2*alpha[2]-theta*w*alpha[3]-w^2*alpha[1]+r+w)*U(xi)+(p^2*alpha[1]+p*q*alpha[3]+q^2*alpha[2])*(diff(diff(U(xi), xi), xi))

Eq, 5.1

> 

ode1111 := collect(ode111/DEtools:-dcoeffs(ode111, U(xi))[1], U(xi))

-U(xi)^3*alpha[4]/(p^2*alpha[1]+p*q*alpha[3]+q^2*alpha[2])+(-theta^2*alpha[2]-theta*w*alpha[3]-w^2*alpha[1]+r+w)*U(xi)/(p^2*alpha[1]+p*q*alpha[3]+q^2*alpha[2])+diff(diff(U(xi), xi), xi)

Note that each term in the Hamiltonian integrand when differentiated wrt U(xi) gives a term in ode. (I have no idea why we do this.) So differentiating wrt xi would give a term in the ode multiplied by diff(U(xi),xi) (using the chain rule). So we can set up a differential equation for the integrand f(xi)

> 

def := diff(f(xi), xi) = ode1111*(diff(U(xi), xi))

diff(f(xi), xi) = (-U(xi)^3*alpha[4]/(p^2*alpha[1]+p*q*alpha[3]+q^2*alpha[2])+(-theta^2*alpha[2]-theta*w*alpha[3]-w^2*alpha[1]+r+w)*U(xi)/(p^2*alpha[1]+p*q*alpha[3]+q^2*alpha[2])+diff(diff(U(xi), xi), xi))*(diff(U(xi), xi))

We want the first integral of this

> 

DEtools:-firint(def, f(xi)); integrand := collect(solve(eval(%, {_C1 = 0, c__1 = 0}), f(xi)), diff, simplify)

f(xi)+((1/4)*alpha[4]*U(xi)^4+(1/2)*(-p^2*alpha[1]-p*q*alpha[3]-q^2*alpha[2])*(diff(U(xi), xi))^2+(1/2)*(theta^2*alpha[2]+theta*w*alpha[3]+w^2*alpha[1]-r-w)*U(xi)^2)/(p^2*alpha[1]+p*q*alpha[3]+q^2*alpha[2])+_C1 = 0

(1/2)*(diff(U(xi), xi))^2-U(xi)^2*(alpha[4]*U(xi)^2+2*theta^2*alpha[2]+2*theta*w*alpha[3]+2*alpha[1]*w^2-2*r-2*w)/(4*p^2*alpha[1]+4*p*q*alpha[3]+4*q^2*alpha[2])

From which we can extract chi and delta, Eq 41.

> 

K := op(1, integrand); P := collect(-op(2, integrand), U(xi))

(1/2)*(diff(U(xi), xi))^2

U(xi)^4*alpha[4]/(4*p^2*alpha[1]+4*p*q*alpha[3]+4*q^2*alpha[2])+(2*theta^2*alpha[2]+2*theta*w*alpha[3]+2*w^2*alpha[1]-2*r-2*w)*U(xi)^2/(4*p^2*alpha[1]+4*p*q*alpha[3]+4*q^2*alpha[2])

Eq. 5.6

> 

H := K+P

(1/2)*(diff(U(xi), xi))^2+U(xi)^4*alpha[4]/(4*p^2*alpha[1]+4*p*q*alpha[3]+4*q^2*alpha[2])+(2*theta^2*alpha[2]+2*theta*w*alpha[3]+2*w^2*alpha[1]-2*r-2*w)*U(xi)^2/(4*p^2*alpha[1]+4*p*q*alpha[3]+4*q^2*alpha[2])

Eq. 5.7

> 

sol1 := U(xi) = m*cos(delta*xi)

U(xi) = m*cos(delta*xi)

Eq. 5.10

> 

eval(convert(H, D), xi = 0); H0 := collect(eval(%, {U(0) = m, (D(U))(0) = 0}), m)

(1/2)*(D(U))(0)^2+U(0)^4*alpha[4]/(4*p^2*alpha[1]+4*p*q*alpha[3]+4*q^2*alpha[2])+(2*theta^2*alpha[2]+2*theta*w*alpha[3]+2*w^2*alpha[1]-2*r-2*w)*U(0)^2/(4*p^2*alpha[1]+4*p*q*alpha[3]+4*q^2*alpha[2])

m^4*alpha[4]/(4*p^2*alpha[1]+4*p*q*alpha[3]+4*q^2*alpha[2])+(2*theta^2*alpha[2]+2*theta*w*alpha[3]+2*w^2*alpha[1]-2*r-2*w)*m^2/(4*p^2*alpha[1]+4*p*q*alpha[3]+4*q^2*alpha[2])

Eq. 5.11 - looks correct

> 

H1 := eval(H, sol1)-H0

(1/2)*m^2*delta^2*sin(delta*xi)^2+m^4*cos(delta*xi)^4*alpha[4]/(4*p^2*alpha[1]+4*p*q*alpha[3]+4*q^2*alpha[2])+(2*theta^2*alpha[2]+2*theta*w*alpha[3]+2*w^2*alpha[1]-2*r-2*w)*m^2*cos(delta*xi)^2/(4*p^2*alpha[1]+4*p*q*alpha[3]+4*q^2*alpha[2])-m^4*alpha[4]/(4*p^2*alpha[1]+4*p*q*alpha[3]+4*q^2*alpha[2])-(2*theta^2*alpha[2]+2*theta*w*alpha[3]+2*w^2*alpha[1]-2*r-2*w)*m^2/(4*p^2*alpha[1]+4*p*q*alpha[3]+4*q^2*alpha[2])

Eq. 5.12 - looks correct

> 

H2 := collect(algsubs(delta*xi = (1/4)*Pi, H1), delta, proc (x) options operator, arrow; collect(x, m) end proc)

(1/4)*m^2*delta^2-(3/16)*m^4*alpha[4]/(p^2*alpha[1]+p*q*alpha[3]+q^2*alpha[2])-(1/16)*(4*theta^2*alpha[2]+4*theta*w*alpha[3]+4*w^2*alpha[1]-4*r-4*w)*m^2/(p^2*alpha[1]+p*q*alpha[3]+q^2*alpha[2])

Eq. 5.13

> 

delta = sort([solve(H2, delta)])[2]; eqphi := `assuming`([simplify(%)], [p^2*alpha[1]+p*q*alpha[3]+q^2*alpha[2] > 0])

delta = (1/2)*((p^2*alpha[1]+p*q*alpha[3]+q^2*alpha[2])*(3*m^2*alpha[4]+4*theta^2*alpha[2]+4*theta*w*alpha[3]+4*w^2*alpha[1]-4*r-4*w))^(1/2)/(p^2*alpha[1]+p*q*alpha[3]+q^2*alpha[2])

delta = (1/2)*(4*alpha[1]*w^2+(4*theta*alpha[3]-4)*w+3*m^2*alpha[4]+4*theta^2*alpha[2]-4*r)^(1/2)/(p^2*alpha[1]+p*q*alpha[3]+q^2*alpha[2])^(1/2)

> 

sol2 := simplify(eval(sol1, eqphi))

U(xi) = m*cos((1/2)*(4*alpha[1]*w^2+(4*theta*alpha[3]-4)*w+3*m^2*alpha[4]+4*theta^2*alpha[2]-4*r)^(1/2)*xi/(p^2*alpha[1]+p*q*alpha[3]+q^2*alpha[2])^(1/2))

> 

simplify(odetest(sol2, ode1111), symbolic)

-4*m*cos((1/2)*(4*alpha[1]*w^2+(4*theta*alpha[3]-4)*w+3*m^2*alpha[4]+4*theta^2*alpha[2]-4*r)^(1/2)*xi/(p^2*alpha[1]+p*q*alpha[3]+q^2*alpha[2])^(1/2))*(cos((1/2)*(4*alpha[1]*w^2+(4*theta*alpha[3]-4)*w+3*m^2*alpha[4]+4*theta^2*alpha[2]-4*r)^(1/2)*xi/(p^2*alpha[1]+p*q*alpha[3]+q^2*alpha[2])^(1/2))^2*m^2*alpha[4]+(3/4)*m^2*alpha[4]+2*alpha[1]*w^2+(2*theta*alpha[3]-2)*w+2*theta^2*alpha[2]-2*r)/(4*p^2*alpha[1]+4*p*q*alpha[3]+4*q^2*alpha[2])

> 

simplify(odetest(sol2, ode111), symbolic)

-(cos((1/2)*(4*alpha[1]*w^2+(4*theta*alpha[3]-4)*w+3*m^2*alpha[4]+4*theta^2*alpha[2]-4*r)^(1/2)*xi/(p^2*alpha[1]+p*q*alpha[3]+q^2*alpha[2])^(1/2))^2*m^2*alpha[4]+(3/4)*m^2*alpha[4]+2*alpha[1]*w^2+(2*theta*alpha[3]-2)*w+2*theta^2*alpha[2]-2*r)*m*cos((1/2)*(4*alpha[1]*w^2+(4*theta*alpha[3]-4)*w+3*m^2*alpha[4]+4*theta^2*alpha[2]-4*r)^(1/2)*xi/(p^2*alpha[1]+p*q*alpha[3]+q^2*alpha[2])^(1/2))

 

 

NULL

Download u-1.mw

seems OK for me...

dharr 9047
March 05 2026
0 0

If I enter www.mapletransaction.org, it loads the page https://mapletransactions.org/index.php/maple and I see the journal home page. This is from an ip address in Canada, but not a university one.

The "about" page says it is a "diamond open access journal" so I don't see why there would be any geo restrictions.

two methods...

dharr 9047
March 02 2026
0 0

@wingho In your listprocedure method, you don't need the psol(0.0001, 0.1); step.

The method I used is the default procedurelist method. The attached does both of them, through to extracting the u(x,t) values. They are both documented on the help page ?pdsolve,numeric

value.mw

workaround...

dharr 9047
February 28 2026
0 0

@nm So a workaround for your version would be to use coulditbe

restart;
assume(_Z1,integer);
coulditbe(-Pi*_Z1+Pi=0);

gives true

remark...

dharr 9047
February 28 2026
0 0

@Alfred_F As I said, I don't have any problem with your solution, but the question is about what odetest does with the ic. You didn't do odetest(sol, [ode, ic], y(x)) - what does that give in Maple 2024?

version?...

dharr 9047
February 27 2026
0 0

@nm So that is different from 2025.2, see worksheet above.

1 2 3 4 5 6 7 Last Page 1 of 98