Advanced Search

dharr

Dr. David Harrington

8205 Reputation

22 Badges

20 years, 337 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 answers submitted by dharr

Bug reporting...

dharr 8205
December 21 2023
0 0

That's quite dramatic! No chugging away for a longtime before losing the kernel as more usually happens. To submit a bug report log in to Mapleprimes, choose the "more" menu, then Submit software change request. In the form there is somewhere to enter additional information, and I often enter the mapleprimes link there (which I copy before I start the form).

Some are there...

dharr 8205
December 17 2023
3 1

The four Appell functions are there, and related functions such as the four Theta functions (JacobiTheta1, etc). For working with q-series and related stuff like Rogers-Ramanujan, see also Frank Garvin's (old-style) Maple packages here. 

> 

restart;
> 

FunctionAdvisor(function_classes);

[trig, trigh, arctrig, arctrigh, elementary, GAMMA_related, Psi_related, Kelvin, Airy, Hankel, Bessel_related, `0F1`, orthogonal_polynomials, Ei_related, erf_related, Kummer, Whittaker, Cylinder, `1F1`, Elliptic_related, Legendre, Chebyshev, `2F1`, Lommel, Struve_related, hypergeometric, Jacobi_related, InverseJacobi_related, Elliptic_doubly_periodic, Weierstrass_related, Zeta_related, complex_components, piecewise_related, Other, Bell, Heun, Appell, trigall, arctrigall, Polylog_related, integral_transforms]

> 

FunctionAdvisor(Appell);

The 4 functions in the "Appell" class are:

[AppellF1, AppellF2, AppellF3, AppellF4]

> 

FunctionAdvisor(definition,AppellF4);

[AppellF4(a, b, c__1, c__2, z__1, z__2) = Sum(Sum(pochhammer(a, _k1+_k2)*pochhammer(b, _k1+_k2)*z__1^_k1*z__2^_k2/(pochhammer(c__1, _k1)*pochhammer(c__2, _k2)*factorial(_k1)*factorial(_k2)), _k2 = 0 .. infinity), _k1 = 0 .. infinity), abs(z__1)^(1/2)+abs(z__2)^(1/2) < 1]

> 

FunctionAdvisor(Jacobi_related);

The 18 functions in the "Jacobi_related" class are:

[JacobiAM, JacobiCD, JacobiCN, JacobiCS, JacobiDC, JacobiDN, JacobiDS, JacobiNC, JacobiND, JacobiNS, JacobiSC, JacobiSD, JacobiSN, JacobiTheta1, JacobiTheta2, JacobiTheta3, JacobiTheta4, JacobiZeta]

> 

FunctionAdvisor(definition,JacobiTheta1);

[JacobiTheta1(z, q) = Sum(2*q^((_k1+1/2)^2)*sin(z*(2*_k1+1))*(-1)^_k1, _k1 = 0 .. infinity), abs(q) < 1]

> 

NULL

Download qFunctions.mw

eta vs Nr...

dharr 8205
December 17 2023
1 0

It's not clear what you want to plot. I'm guessing eta vs Nr for Nr=0.1, 0.3, 0.6, (for x=0) so there are only 3 points on this plot?

FE-2.mw

 

Inconsistent system...

dharr 8205
December 17 2023
3 1

Since your equations are linear, I converted them to Matrix form, then LinearSolve reports that the system is inconsistent. Note that the determinant of your matrix is zero and the rank is only 12. (I assume that eq[16] being the only one with a constant term is correct - but is the constant correct?) So there are relationships between the equations, and perhaps between the original differential equations. Using Nullspace to see what the relationships are between columns of the Matrix doesn't show anything obvious.

Perhaps Maple can solve the system of equations symbolically so it is easier to see what the relationships might be. Or you could set up the differential equations in Matrix form to make it more obvious what is happening.

Timoshenko_Beam.mw

another way for display only...

dharr 8205
December 16 2023
2 0

This prints in the same typeface as usual.

> 

restart

Setup so pp(expr) prettyprints expr with "," as decimal separator

> 

pp:=(expr)->Typesetting:-Typeset(expr):
`print/mn`:=x->Typesetting:-mn(StringTools:-Subs("."=",",x)):
> 

q := 1.23+sqrt((1.5+y)/sqrt(sin(1.45)+3.5*x))
pp(%)

Parse:-ConvertTo1D, "invalid input %1", Typesetting:-mprintslash([Typesetting:-mrow(Typesetting:-mn("1,23"), Typesetting:-mo("&plus;"), Typesetting:-msqrt(Typesetting:-mfrac(Typesetting:-mrow(Typesetting:-mn("1,5"), Typesetting:-mo("&plus;"), Typesetting:-mi("y")), Typesetting:-msqrt(Typesetting:-mrow(Typesetting:-mn("0,9927129910"), Typesetting:-mo("&plus;"), Typesetting:-mrow(Typesetting:-mn("3,5"), Typesetting:-mo("&InvisibleTimes;"), Typesetting:-mi("x")))))))], [Typesetting:-mrow(Typesetting:-mn("1.23"), Typesetting:-mo("&plus;"), Typesetting:-msqrt(Typesetting:-mfrac(Typesetting:-mrow(Typesetting:-mn("1.5"), Typesetting:-mo("&plus;"), Typesetting:-mi("y")), Typesetting:-msqrt(Typesetting:-mrow(Typesetting:-mn("0.9927129910"), Typesetting:-mo("&plus;"), Typesetting:-mrow(Typesetting:-mn("3.5"), Typesetting:-mo("&InvisibleTimes;"), Typesetting:-mi("x")))))))])

> 

NULL

Download comma.mw

evaluation issue...

dharr 8205
December 16 2023
3 0

I fixed this in a quick-and-dirty way by adding eval at some places (shown in red). If the logic flow changes, you will need more, but it is better to adjust the way you are progamming. You assign expressions with xi in it to variables and then do xi:=1 or similar. Within a procedure there is only one-level evaluation so this does not automatically change the value of the variables. You should use say, xi1:=1, then eval(something with xi in it, xi = xi1).

SiS_Quadratic_(goal).mw 

(I only tested to the line FirmModelH(0.8, 0.15, 0.15, 0.10)[3])

Check how this simpler code works:

test:=proc(a) local xi, q, q2;
  q:=a*xi;
  xi:=2;
  q2:=eval(q);
  return q,q2;
end proc:

test(3);

gives 

plot...

dharr 8205
December 11 2023
0 0

Not sure I entirely understand the details of what you want. I used blue axes, but if you want arrow axes with better-positioned labels you could use axes=none and make your own with arrow and textplot3d. (Proportions look better in the worksheet.)

``

> 

restart
> 

with(plots); with(plottools)
> 

``

> 

cuboid1 := cuboid([0, -2, -4], [10, 2, 4], color = "red", style = wireframe); ell := textplot3d([5, 0, -4.5, "&ell;"], align = below, color = "blue", font = [Times, italic, 20]); text_inside := textplot3d([5, 0, 0, typeset(Q(x, z, t) = Q[0]*delta(-t*v+x)*sin(p*z))], align = below, color = "green", font = [Times, italic, 16]); display(cuboid1, text_inside, ell, scaling = constrained, axes = normal, axis = [color = "blue"], view = [-1 .. 12, -6 .. 6, -6 .. 6], orientation = [-80, 80], labels = [x, y, z], tickmarks = [0, 0, 0], labelfont = [Times, italic, 16])

 

NULL

Download rectangular_cylinder.mw

solve whole system...

dharr 8205
December 11 2023
1 3

The whole system of all 8 equations can all be solved at once. You did not have a value for the "omegat" parameter, so I arbitrarily set that to 1. You also had a error in setting values for parameters at the beginning of the sheet, but since that line was redundant, I just removed it.

 New.mw

parity of rotations...

dharr 8205
December 09 2023
2 0

I didn't track it down exactly, but it is related the the fact that rotations of [1,2,...,n] considered as permutations all have even parity if n is odd, but have alternating odd, even, parity if n is even (see in attached). The transpose of the adjoint is the matrix of cofactors, which are (-1)^(i+j) times the determinant of the minor, and the signs of the terms in the determinant have signs depending on the parity of the permutation of the indices. So in using rotations sometimes you have alternating signs and don't need your matrix P for the (-1)^(i+j) , and sometimes you don't have alternating signs, so need P. 

The attached is mainly me figuring out what you were doing by translating to the language of permutation matrices.

Adjoint.mw

Using plot component...

dharr 8205
December 06 2023
3 9

I think @acer mentioned using a plot component. Here is an implementation of that. Run the worksheet then follow the instructions. There is code in the reset button click event handler and the plot click event handler.

Rename DealeyPlaza.txt to DealeyPlaza.jpg

DealeyPlaza.txt

MeasureImage.mw

integer solutions using algcurves...

dharr 8205
December 06 2023
4 0

Find the smallest positive integer solution to eq below.

> 

restart
> 

with(algcurves)
> 

eq := a/(b+c)+b/(a+c)+c/(a+b) = 4

a/(b+c)+b/(a+c)+c/(a+b) = 4

Smallest (!) positive solution is

> 

soln1 := {a = 4373612677928697257861252602371390152816537558161613618621437993378423467772036, b = 36875131794129999827197811565225474825492979968971970996283137471637224634055579, c = 154476802108746166441951315019919837485664325669565431700026634898253202035277999}; eval(eq, soln1)

4 = 4

We more-or-less follow the method given in https://www.agftutoring.com/x-yz-y-xz-z-xy-4/
Rearrange as a polynomial

> 

normal(eq); eq2 := expand(numer(lhs(%))-rhs(%)*denom(lhs(%)))

a^3-3*a^2*b-3*a^2*c-3*a*b^2-5*a*b*c-3*a*c^2+b^3-3*b^2*c-3*b*c^2+c^3

Convert this homogeneous polynomial to a bivariate polynomial - elliptic curves have genus = 1

> 

f := eval(eq2, c = 1); genus(f, a, b)

a^3-3*a^2*b-3*a*b^2+b^3-3*a^2-5*a*b-3*b^2-3*a-3*b+1

1

The j_invariant has the same value despite the various transformations (birational equivalence)

> 

j_invariant(f, a, b)

1408317602329/2153060

Maple finds a Weierstrass form

> 

wform := Weierstrassform(f, a, b, X, Y); w1 := wform[1]

X^3-(11209/48)*X-1185157/864+Y^2

That is different from the one in the reference, but that one may be transformed to w1

> 

w3 := y^2-x^3-109*x^2-224*x; j_invariant(w3, x, y); simplify(eval((1/64)*w3, {x = -4*X-109/3, y = 8*Y}))

-x^3-109*x^2+y^2-224*x

1408317602329/2153060

X^3-(11209/48)*X-1185157/864+Y^2

However, for the addition routine below we would like to have it in the standard form y^2 = x^3+A*x+B, with A and B integers.

Try simply scaling X and Y and the whole expression; we also need to change the sign of X.

> 

trans := {X = -x/6^2, Y = y/(72*3)}; w := eval(864*(9*6)*w1, trans)

{X = -(1/36)*x, Y = (1/216)*y}

-x^3+y^2+302643*x-63998478

> 

j_invariant(w, x, y)

1408317602329/2153060

Conversion between f and w

> 

itrans := solve(trans, {x, y}); ab_to_xy := eval(itrans, {X = wform[2], Y = wform[3]})

{x = -36*X, y = 216*Y}

{x = 216*a^2-864*a*b+216*b^2-612*a-612*b-465, y = 7776*a^3-31104*a^2*b+7776*a*b^2-22680*a^2-19440*a*b-648*b^2-25164*a+1836*b+3888}

Test with a known solution (can be found using method below)

> 

xy := eval(ab_to_xy, {a = 4/11, b = -1/11}); eval(w, %)

{x = -573, y = -7020}

0

Conversion between w and f

> 

xy_to_ab := simplify(eval({a = wform[4], b = wform[5]}, trans)); `union`(eval(xy_to_ab, xy), {c = 1}); eval(eq, %)

{a = (3*x-2493+y)/(-10260+36*x), b = (3*x-2493-y)/(-10260+36*x)}

{a = 4/11, b = -1/11, c = 1}

4 = 4

Find the discriminant of the cubic

> 

Cubic := rhs(isolate(w, y^2)); d := discrim(Cubic, x); ifactor(d)

x^3-302643*x+63998478

292921436021760

``(2)^10*``(3)^12*``(5)*``(7)^2*``(13)^3

To find the first solution, we want either y = 0, or a y value for which y^2 divides the disciminant, for which the corresponding x is an integer  (Nagell-Lutz) .
So we test all the divisors of 7*2^5*3^6

> 

divs := `union`({0}, NumberTheory:-Divisors(7*2^5*3^6))
> 

j := 0; for i in divs do facs := factors(eval(w, y = i))[2]; if nops(facs) > 1 then j := j+1; sol[j] := {x = solve(select(proc (z) options operator, arrow; type(z[1], linear) end proc, facs)[1][1], x), y = i}; print(sol[j]) end if end do

{x = 327, y = 0}

{x = 291, y = 756}

We found two solutions - convert back to a,b,c form and check. First one is OK for eq2, but not for the original eq because a denominator is zero

> 

`union`(eval(xy_to_ab, sol[1]), {c = 1}); eval(eq2, %); eval(eq, `%%`)

{a = -1, b = -1, c = 1}

0

Error, numeric exception: division by zero

So the only solution is the second one. (The torsion subgroup has only one generator.)

> 

`union`(eval(xy_to_ab, sol[2]), {c = 1}); eval(eq2, %); eval(eq, `%%`)

{a = -4, b = -11, c = 1}

0

4 = 4

To find other points we combine existing point(s) until we find a positive one
Define addition of points in the group algebra - see Wikipedia

A is the coefficient in y^2=x^3+A*x+B

> 

algadd := proc(P::[rational,rational], Q::[rational,rational], A::rational)
  local s,xP,yP,xQ,yQ,xR,yR;
  xP,yP := P[];
  xQ,yQ := Q[];
  if xP <> xQ then
    s := (yP - yQ)/(xP - xQ);
    xR := s^2 - xP - xQ;
    yR := yP - s*(xP - xR);
    [xR, -yR];
  elif yP = -yQ then
    [0, 0] # infinity point not representable
  else # P = Q,sum is tangent
    s := (3*xP^2 + A)/(2*yP);
    xR := s^2 - 2*xP;
    yR := yP - s*(xP - xR);
    [xR, -yR]
  end if;
end proc:
> 

A := coeff(Cubic, x, 1)

-302643

> 

P[1] := eval([x, y], sol[2])

[291, 756]

Incrementally add P[1] until we get a positive solution

> 

for i to 8 do P[i+1] := algadd(P[i], P[1], A); soln := eval(xy_to_ab, `~`[`=`]([x, y], P[i+1])) end do
> 

soln2 := {a = numer(eval(a, soln)), b = numer(eval(b, soln)), c = denom(eval(a, soln))}

{a = 36875131794129999827197811565225474825492979968971970996283137471637224634055579, b = 154476802108746166441951315019919837485664325669565431700026634898253202035277999, c = 4373612677928697257861252602371390152816537558161613618621437993378423467772036}

> 

eval(eq, soln2)

4 = 4

> 

NULL

Download elliptic.mw

dcoeffs...

dharr 8205
December 04 2023
2 14

dcoeffs can do this - the order is determined by the type of derivative, not the way the ode is written; they are the same in your example so hope this generalizes the way you want.

NULL

> 

restart
> 

with(PDEtools)
> 

ode := diff(y(x), x, x) = c*(diff(y(x), x))-2*a*y(x)+b*x

diff(diff(y(x), x), x) = c*(diff(y(x), x))-2*a*y(x)+b*x

> 

cfs := [dcoeffs(rhs(ode)-lhs(ode), y(x))]; vars := [seq(cat(A, i), i = 1 .. 4)]; add(`~`[`*`](cfs, vars)); poly := subs(x = 1, %)

[-1, c, -2*a, b*x]

[A1, A2, A3, A4]

A4*b*x+A2*c-2*A3*a-A1

A2*c-2*A3*a+A4*b-A1

> 

NULL

NULL

Download dcoeffs.mw

tiny values...

dharr 8205
November 30 2023
1 0
> 

restart
> 

with(Student[MultivariateCalculus]):
> 

f := 39584968580329795728950940517214770307434335*xx^3*yy^6*(1/21778071482940061661655974875633165533184)+2207379816207475241162406248223006569040862935*xx^5*yy^8*(1/2787593149816327892691964784081045188247552)+47950825635610780986659544491454706340397108297*xx^5*yy^7*(1/22300745198530623141535718272648361505980416)-588774433706353379897742534304221654039246663*xx^5*yy^6*(1/348449143727040986586495598010130648530944)-44608078263668464626393951292252447406629869273*xx^5*yy^5*(1/22300745198530623141535718272648361505980416)+1613038118657167505912389296857854524947676825*xx^5*yy^4*(1/1393796574908163946345982392040522594123776)+23570688854853763073042723518782612790921757535*xx^5*yy^3*(1/22300745198530623141535718272648361505980416)-113510140727511300460098712979462156361337425*xx^5*yy^2*(1/348449143727040986586495598010130648530944)-6817973449093402642853212701104432585928821163*xx^5*yy*(1/22300745198530623141535718272648361505980416)+184838927094446995029201369223921105703104647*xx^5*(1/2787593149816327892691964784081045188247552)+31380186488931551370058361496245928395816772575*xx^4*yy^9*(1/1393796574908163946345982392040522594123776)-1758702445038817232726176779731884586549332868025*xx^4*yy^8*(1/44601490397061246283071436545296723011960832)-33218490572036542393092937176469859040906121155*xx^4*yy^7*(1/348449143727040986586495598010130648530944)+116540829629507365267125159526451609264014215*xx^7*yy^6*(1/87112285931760246646623899502532662132736)-2146702909675882809503682033933399905335826325*xx^9*yy^9*(1/11150372599265311570767859136324180752990208)+1587967252519403636411870604735180043125989625*xx^9*yy^8*(1/5575186299632655785383929568162090376495104)-20361225581568567923686744589522827658576624955*xx^3*yy^7*(1/11150372599265311570767859136324180752990208)-11540959773500599403794316292492996114189538863*xx^2*yy^2*(1/5575186299632655785383929568162090376495104)-1174244552874873223035231031480900497934023075*xx^3*yy^8*(1/1393796574908163946345982392040522594123776)+27287439738914744607616926917914225474665410565*xx^2*yy^3*(1/174224571863520493293247799005065324265472)-206512033439850904054937113093163624192322042825*xx^8*yy^4*(1/44601490397061246283071436545296723011960832)-62755544772437504320590342390381422715234113715/89202980794122492566142873090593446023921664+11170081785792631086653879206603595320491089331*xx^7*yy^5*(1/11150372599265311570767859136324180752990208)-9205355621994819342146712860571987786619361601*xx*yy*(1/44601490397061246283071436545296723011960832)-652299342907430898149182084981866414949696905*xx^7*yy^4*(1/696898287454081973172991196020261297061888)+27484692689867334306687311759874973819976026005*xx*yy^3*(1/44601490397061246283071436545296723011960832)-4303517165264733669855129139552505045324631645*xx^7*yy^3*(1/11150372599265311570767859136324180752990208)+10578825782023300845453772557509072093336001*xx^7*yy^2*(1/43556142965880123323311949751266331066368)+157001869330425518481531763580902779395436599415*xx^8*yy^6*(1/22300745198530623141535718272648361505980416)-930314746723434588666177195703059675161177190255*xx^2*yy^6*(1/5575186299632655785383929568162090376495104)-3917684154726736823398471536296978037714283086195*yy^8*(1/89202980794122492566142873090593446023921664)+998213736763384913910074759047227544847506773*xx^7*yy*(1/11150372599265311570767859136324180752990208)+1583056855557692418384969876461998197073089695*xx*yy^4*(1/2787593149816327892691964784081045188247552)-104255809907916433055923335622932126645726549*xx*yy^2*(1/696898287454081973172991196020261297061888)+1970986683407627074325019523003479974617451789943*yy^6*(1/22300745198530623141535718272648361505980416)-77131555128675321096947207038878222843991869993*yy^7*(1/696898287454081973172991196020261297061888)-29946355461657315300256240552185966952551471*xx^7*(1/1393796574908163946345982392040522594123776)-125283292999146417157156696376640452081866835*xx^3*(1/1393796574908163946345982392040522594123776)+23458516464006675395891679247259419002768896835*xx^3*yy^5*(1/11150372599265311570767859136324180752990208)-539977758872163289054492124375185771143918033*xx^6*yy*(1/696898287454081973172991196020261297061888)-3479476522267890993628796487849129439635143625*xx^8*yy^3*(1/696898287454081973172991196020261297061888)+18712604797880071317805036942199122521197359575*xx^8*yy^2*(1/22300745198530623141535718272648361505980416)+43423414494451507811145033075147441881593811799*yy^2*(1/22300745198530623141535718272648361505980416)-15637727799880882327290754576104647826715168925*xx^3*yy^3*(1/11150372599265311570767859136324180752990208)+1206817075246069632318716986669541278160772775*xx^9*yy^4*(1/2787593149816327892691964784081045188247552)-2038600361316622246653155899145012259420048867785*yy^4*(1/44601490397061246283071436545296723011960832)+30423874459994412977383604476886160940746185*xx^9*(1/5575186299632655785383929568162090376495104)+220816865194317615868568855814620996552449073*xx*(1/5575186299632655785383929568162090376495104)+5138909461003175489938484170634052266819688725*xx^9*yy^3*(1/44601490397061246283071436545296723011960832)+29341459645317546529685572705520876577051855*xx^3*yy^2*(1/87112285931760246646623899502532662132736)-1277356081222180962342283013232991241852904465*xx^6*yy^9*(1/696898287454081973172991196020261297061888)+17639360745426635511855086638766468926126459875*xx^9*yy^7*(1/44601490397061246283071436545296723011960832)+15350689937843699961175740256400109996121380375*xx^8*yy^5*(1/1393796574908163946345982392040522594123776)-285743684916570536194588196441080828723328178675*xx^8*yy^8*(1/89202980794122492566142873090593446023921664)+1869246621670048362557342074310025153518449965*xx^8*yy*(1/2787593149816327892691964784081045188247552)-2255097230860381206152749351617455809672044745*xx*yy^9*(1/11150372599265311570767859136324180752990208)-4089215965643055747590786827106386135115380275*xx^8*(1/89202980794122492566142873090593446023921664)+35122173917479363738100862234581108137514304171*xx^2*(1/22300745198530623141535718272648361505980416)-72716798311978341010558827315982986191821905*xx^9*yy^2*(1/696898287454081973172991196020261297061888)+57447439083834576362467553225131370438848237035*xx^6*yy^8*(1/22300745198530623141535718272648361505980416)-431284328058774504067793959976795724976545555*xx^9*yy^6*(1/696898287454081973172991196020261297061888)-17449701902039745490242163912540688306429882361*xx^2*yy*(1/696898287454081973172991196020261297061888)-56566850002827011453690682806041619180254985625*yy^3*(1/696898287454081973172991196020261297061888)-1197236208181378637639504269592639035279087665*xx^9*yy*(1/44601490397061246283071436545296723011960832)+216255546256559295251079313253452049445763455*xx^7*yy^9*(1/348449143727040986586495598010130648530944)+76828297887427851822683521168415270943435162685*yy^9*(1/2787593149816327892691964784081045188247552)+36390552938954376406834468187448925576623439893*xx^2*yy^7*(1/174224571863520493293247799005065324265472)+8094790880015327525694605814920739418439287725*xx^8*yy^9*(1/2787593149816327892691964784081045188247552)+941109349474535911451616661821106567867537125*xx^3*yy^9*(1/1393796574908163946345982392040522594123776)-5023626067733175609651265492842895195168362165*xx^5*yy^9*(1/5575186299632655785383929568162090376495104)-36304948749180317956941914133403396762716230691*xx*yy^5*(1/44601490397061246283071436545296723011960832)-12993287722661922638788467553649639108437064835*xx^9*yy^5*(1/44601490397061246283071436545296723011960832)+211134394987302797546644924545169826774270265159*yy^5*(1/1393796574908163946345982392040522594123776)+35696532930567486560276536615522532283474689213*yy*(1/2787593149816327892691964784081045188247552)-48412290717709997717153300332089796247538326265*xx^4*(1/44601490397061246283071436545296723011960832)+1173296429365947392287371443632107462978009165*xx^6*yy^7*(1/174224571863520493293247799005065324265472)+17196469545705046799299985950707233685621881055*xx^4*yy*(1/1393796574908163946345982392040522594123776)+929769947314964740179937673332890647768037984465*xx^2*yy^4*(1/11150372599265311570767859136324180752990208)-868641325364973493898126340263842300348545855*xx^7*yy^8*(1/1393796574908163946345982392040522594123776)+5011420945327438626354964312196465908094234685*xx^3*yy*(1/11150372599265311570767859136324180752990208)-9551461763890264957289963973620923748598225435*xx^4*yy^2*(1/11150372599265311570767859136324180752990208)-23673134207774883972271882396704370580007933039*xx^6*yy^6*(1/5575186299632655785383929568162090376495104)+5795161625895678368156852916105373987594511979*xx^6*(1/22300745198530623141535718272648361505980416)-2937701213452088192123555543440803264914467299*xx^6*yy^5*(1/348449143727040986586495598010130648530944)-14785537121406447202257499440081382142298519099*xx^7*yy^7*(1/11150372599265311570767859136324180752990208)+14159347676475748959036290080103848146860867025*xx^6*yy^4*(1/11150372599265311570767859136324180752990208)-6686861200533386632065997818427854246215113305*xx^8*yy^7*(1/696898287454081973172991196020261297061888)-26051472095770585704126329008135447818638784275*xx^4*yy^3*(1/348449143727040986586495598010130648530944)+749877940244270735637721966049124917356845885*xx^6*yy^3*(1/174224571863520493293247799005065324265472)+27046038795224386955728969793334632924015008227*xx*yy^7*(1/44601490397061246283071436545296723011960832)+782685832362921584689673760969891945953777553*xx^6*yy^2*(1/5575186299632655785383929568162090376495104)-100809382380090436397261413740272360141145204891*xx^2*yy^5*(1/348449143727040986586495598010130648530944)+2168816628024980374461014350770096009019357665*xx*yy^8*(1/5575186299632655785383929568162090376495104)-765302392604646459013613426858243443467023490875*xx^4*yy^4*(1/22300745198530623141535718272648361505980416)+1872760743346397986120124413411813119412045269675*xx^2*yy^8*(1/22300745198530623141535718272648361505980416)+94251624724512021502035994822030873708141367565*xx^4*yy^5*(1/696898287454081973172991196020261297061888)-35643509355104072817665294345590475660747146425*xx^2*yy^9*(1/696898287454081973172991196020261297061888)+843981485493394825713526892530506348990296828805*xx^4*yy^6*(1/11150372599265311570767859136324180752990208)-590212436135125327923049635849260481403670583*xx*yy^6*(1/696898287454081973172991196020261297061888)-851688199122087410134053760306093104684621525*xx^3*yy^4*(1/696898287454081973172991196020261297061888)
> 

g := .5*(1+tanh(f))

Values are close to zero, so need to scale function to reasonable values

> 

TaylorApproximation(evalf(g), [xx, yy] = [0, 0], 6)

0.2259369109e-600*xx^6+0.1496133470e-597*xx^5*yy+0.4645575308e-595*xx^4*yy^2+0.7998759715e-593*xx^3*yy^3+0.9696415420e-591*xx^2*yy^4+0.6264381195e-589*xx*yy^5+0.3379234176e-587*yy^6+0.5849644085e-602*xx^5+0.3631645737e-599*xx^4*yy+0.9376535155e-597*xx^3*yy^2+0.1515139614e-594*xx^2*yy^3+0.1223177462e-592*xx*yy^4+0.7916210715e-591*yy^5+0.1419522894e-603*xx^4+0.7327837035e-601*xx^3*yy+0.1775662734e-598*xx^2*yy^2+0.1910714715e-596*xx*yy^3+0.1545397580e-594*yy^4+0.2863414606e-605*xx^3+0.1387337376e-602*xx^2*yy+0.2238559204e-600*xx*yy^2+0.2413560442e-598*yy^3+0.5419743765e-607*xx^2+0.1748457882e-604*xx*yy+0.2827108662e-602*yy^2+0.6828367775e-609*xx+0.2207703247e-606*yy

> 

scale := 10^(-610)
> 

h := evalf(TaylorApproximation(g/scale, [xx, yy] = [0, 0], 35))
> 

plot3d(h, xx = -1 .. 1, yy = -10 .. 10, color = red, style = surface)

> 

 

NULL

Download taylorProblem.mw

LinearSolve...

dharr 8205
November 30 2023
2 0

This can be done with LinearSolve as shown. c1 etc are very complicated expressions that perhaps can be simplified, but I did not wait to see if that were the case. This method assumes that the equations are independent (Determinant(A)<>0), and pays no attention to any relationships that might exist between the Bessel functions.

NULL

> 

restart;
> 

with(LinearAlgebra):
> 

fn1 := (2*BesselK(1, alpha1*sigma)*T1*alpha1^3*b1*c+2*BesselK(1, alpha1*sigma)*T1*T2*alpha1*b1*c+BesselK(1, alpha1*sigma)*alpha1*b1*c+BesselK(1, alpha1*sigma)*alpha1*b1+BesselK(0, alpha1*sigma))*c1+(-2*BesselI(1, alpha1*sigma)*T1*alpha1^3*b1*c-2*BesselI(1, alpha1*sigma)*T1*T2*alpha1*b1*c-BesselI(1, alpha1*sigma)*alpha1*b1*c-BesselI(1, alpha1*sigma)*alpha1*b1+BesselI(0, alpha1*sigma))*c2+(2*BesselK(1, alpha2*sigma)*T1*alpha2^3*b1*c+2*BesselK(1, alpha2*sigma)*T1*T2*alpha2*b1*c+BesselK(1, alpha2*sigma)*alpha2*b1*c+BesselK(1, alpha2*sigma)*alpha2*b1+BesselK(0, alpha2*sigma))*c3+(-2*BesselI(1, alpha2*sigma)*T1*alpha2^3*b1*c-2*BesselI(1, alpha2*sigma)*T1*T2*alpha2*b1*c-BesselI(1, alpha2*sigma)*alpha2*b1*c-BesselI(1, alpha2*sigma)*alpha2*b1+BesselI(0, alpha2*sigma))*c4+ur/alpha^2+(-2*A1*K*BesselI(1, k*sigma)*T1*b1*c*k^3+2*B1*K*BesselK(1, k*sigma)*T1*b1*c*k^3-2*A1*K*BesselI(1, k*sigma)*T1*T2*b1*c*k+2*B1*K*BesselK(1, k*sigma)*T1*T2*b1*c*k-2*A1*BesselI(1, k*sigma)*T1*b1*c*k*k1+2*B1*BesselK(1, k*sigma)*T1*b1*c*k*k1-A1*K*BesselI(1, k*sigma)*b1*c*k+B1*K*BesselK(1, k*sigma)*b1*c*k-A1*K*BesselI(1, k*sigma)*b1*k+B1*K*BesselK(1, k*sigma)*b1*k+A1*K*BesselI(0, k*sigma)+B1*K*BesselK(0, k*sigma))*EZ;

(2*BesselK(1, alpha1*sigma)*T1*alpha1^3*b1*c+2*BesselK(1, alpha1*sigma)*T1*T2*alpha1*b1*c+BesselK(1, alpha1*sigma)*alpha1*b1*c+BesselK(1, alpha1*sigma)*alpha1*b1+BesselK(0, alpha1*sigma))*c1+(-2*BesselI(1, alpha1*sigma)*T1*alpha1^3*b1*c-2*BesselI(1, alpha1*sigma)*T1*T2*alpha1*b1*c-BesselI(1, alpha1*sigma)*alpha1*b1*c-BesselI(1, alpha1*sigma)*alpha1*b1+BesselI(0, alpha1*sigma))*c2+(2*BesselK(1, alpha2*sigma)*T1*alpha2^3*b1*c+2*BesselK(1, alpha2*sigma)*T1*T2*alpha2*b1*c+BesselK(1, alpha2*sigma)*alpha2*b1*c+BesselK(1, alpha2*sigma)*alpha2*b1+BesselK(0, alpha2*sigma))*c3+(-2*BesselI(1, alpha2*sigma)*T1*alpha2^3*b1*c-2*BesselI(1, alpha2*sigma)*T1*T2*alpha2*b1*c-BesselI(1, alpha2*sigma)*alpha2*b1*c-BesselI(1, alpha2*sigma)*alpha2*b1+BesselI(0, alpha2*sigma))*c4+ur/alpha^2+(-2*A1*K*BesselI(1, k*sigma)*T1*b1*c*k^3+2*B1*K*BesselK(1, k*sigma)*T1*b1*c*k^3-2*A1*K*BesselI(1, k*sigma)*T1*T2*b1*c*k+2*B1*K*BesselK(1, k*sigma)*T1*T2*b1*c*k-2*A1*BesselI(1, k*sigma)*T1*b1*c*k*k1+2*B1*BesselK(1, k*sigma)*T1*b1*c*k*k1-A1*K*BesselI(1, k*sigma)*b1*c*k+B1*K*BesselK(1, k*sigma)*b1*c*k-A1*K*BesselI(1, k*sigma)*b1*k+B1*K*BesselK(1, k*sigma)*b1*k+A1*K*BesselI(0, k*sigma)+B1*K*BesselK(0, k*sigma))*EZ

> 

fn2 := (-2*BesselK(1, alpha1)*T1*alpha1^3*b1*c-2*BesselK(1, alpha1)*T1*T2*alpha1*b1*c-BesselK(1, alpha1)*alpha1*b1*c-BesselK(1, alpha1)*alpha1*b1+BesselK(0, alpha1))*c1+(2*BesselI(1, alpha1)*T1*alpha1^3*b1*c+2*BesselI(1, alpha1)*T1*T2*alpha1*b1*c+BesselI(1, alpha1)*alpha1*b1*c+BesselI(1, alpha1)*alpha1*b1+BesselI(0, alpha1))*c2+(-2*BesselK(1, alpha2)*T1*alpha2^3*b1*c-2*BesselK(1, alpha2)*T1*T2*alpha2*b1*c-BesselK(1, alpha2)*alpha2*b1*c-BesselK(1, alpha2)*alpha2*b1+BesselK(0, alpha2))*c3+(2*BesselI(1, alpha2)*T1*alpha2^3*b1*c+2*BesselI(1, alpha2)*T1*T2*alpha2*b1*c+BesselI(1, alpha2)*alpha2*b1*c+BesselI(1, alpha2)*alpha2*b1+BesselI(0, alpha2))*c4+ur/alpha^2+(2*A1*K*BesselI(1, k)*T1*b1*c*k^3-2*B1*K*BesselK(1, k)*T1*b1*c*k^3+2*A1*K*BesselI(1, k)*T1*T2*b1*c*k-2*B1*K*BesselK(1, k)*T1*T2*b1*c*k+2*A1*BesselI(1, k)*T1*b1*c*k*k1-2*B1*BesselK(1, k)*T1*b1*c*k*k1+A1*K*BesselI(1, k)*b1*c*k-B1*K*BesselK(1, k)*b1*c*k+A1*K*BesselI(1, k)*b1*k-B1*K*BesselK(1, k)*b1*k+A1*K*BesselI(0, k)+B1*K*BesselK(0, k))*EZ;

(-2*BesselK(1, alpha1)*T1*alpha1^3*b1*c-2*BesselK(1, alpha1)*T1*T2*alpha1*b1*c-BesselK(1, alpha1)*alpha1*b1*c-BesselK(1, alpha1)*alpha1*b1+BesselK(0, alpha1))*c1+(2*BesselI(1, alpha1)*T1*alpha1^3*b1*c+2*BesselI(1, alpha1)*T1*T2*alpha1*b1*c+BesselI(1, alpha1)*alpha1*b1*c+BesselI(1, alpha1)*alpha1*b1+BesselI(0, alpha1))*c2+(-2*BesselK(1, alpha2)*T1*alpha2^3*b1*c-2*BesselK(1, alpha2)*T1*T2*alpha2*b1*c-BesselK(1, alpha2)*alpha2*b1*c-BesselK(1, alpha2)*alpha2*b1+BesselK(0, alpha2))*c3+(2*BesselI(1, alpha2)*T1*alpha2^3*b1*c+2*BesselI(1, alpha2)*T1*T2*alpha2*b1*c+BesselI(1, alpha2)*alpha2*b1*c+BesselI(1, alpha2)*alpha2*b1+BesselI(0, alpha2))*c4+ur/alpha^2+(2*A1*K*BesselI(1, k)*T1*b1*c*k^3-2*B1*K*BesselK(1, k)*T1*b1*c*k^3+2*A1*K*BesselI(1, k)*T1*T2*b1*c*k-2*B1*K*BesselK(1, k)*T1*T2*b1*c*k+2*A1*BesselI(1, k)*T1*b1*c*k*k1-2*B1*BesselK(1, k)*T1*b1*c*k*k1+A1*K*BesselI(1, k)*b1*c*k-B1*K*BesselK(1, k)*b1*c*k+A1*K*BesselI(1, k)*b1*k-B1*K*BesselK(1, k)*b1*k+A1*K*BesselI(0, k)+B1*K*BesselK(0, k))*EZ

> 

fn3 := -T1*alpha1*(alpha1^2+T2)*(BesselK(0, alpha1*sigma)*alpha1*b2*sigma+BesselK(1, alpha1*sigma)*b2*s1+BesselK(1, alpha1*sigma)*b2+BesselK(1, alpha1*sigma)*sigma)*c1/sigma+T1*alpha1*(alpha1^2+T2)*(-BesselI(0, alpha1*sigma)*alpha1*b2*sigma+BesselI(1, alpha1*sigma)*b2*s1+BesselI(1, alpha1*sigma)*b2+BesselI(1, alpha1*sigma)*sigma)*c2/sigma-T1*alpha2*(alpha2^2+T2)*(BesselK(0, alpha2*sigma)*alpha2*b2*sigma+BesselK(1, alpha2*sigma)*b2*s1+BesselK(1, alpha2*sigma)*b2+BesselK(1, alpha2*sigma)*sigma)*c3/sigma+T1*alpha2*(alpha2^2+T2)*(-BesselI(0, alpha2*sigma)*alpha2*b2*sigma+BesselI(1, alpha2*sigma)*b2*s1+BesselI(1, alpha2*sigma)*b2+BesselI(1, alpha2*sigma)*sigma)*c4/sigma+T1*EZ*k*(K*k^2+K*T2+k1)*(-A1*BesselI(0, k*sigma)*b2*k*sigma-B1*BesselK(0, k*sigma)*b2*k*sigma+A1*BesselI(1, k*sigma)*b2*s1-B1*BesselK(1, k*sigma)*b2*s1+A1*BesselI(1, k*sigma)*b2+A1*BesselI(1, k*sigma)*sigma-B1*BesselK(1, k*sigma)*b2-B1*BesselK(1, k*sigma)*sigma)/sigma;

-T1*alpha1*(alpha1^2+T2)*(BesselK(0, alpha1*sigma)*alpha1*b2*sigma+BesselK(1, alpha1*sigma)*b2*s1+BesselK(1, alpha1*sigma)*b2+BesselK(1, alpha1*sigma)*sigma)*c1/sigma+T1*alpha1*(alpha1^2+T2)*(-BesselI(0, alpha1*sigma)*alpha1*b2*sigma+BesselI(1, alpha1*sigma)*b2*s1+BesselI(1, alpha1*sigma)*b2+BesselI(1, alpha1*sigma)*sigma)*c2/sigma-T1*alpha2*(alpha2^2+T2)*(BesselK(0, alpha2*sigma)*alpha2*b2*sigma+BesselK(1, alpha2*sigma)*b2*s1+BesselK(1, alpha2*sigma)*b2+BesselK(1, alpha2*sigma)*sigma)*c3/sigma+T1*alpha2*(alpha2^2+T2)*(-BesselI(0, alpha2*sigma)*alpha2*b2*sigma+BesselI(1, alpha2*sigma)*b2*s1+BesselI(1, alpha2*sigma)*b2+BesselI(1, alpha2*sigma)*sigma)*c4/sigma+T1*EZ*k*(K*k^2+K*T2+k1)*(-A1*BesselI(0, k*sigma)*b2*k*sigma-B1*BesselK(0, k*sigma)*b2*k*sigma+A1*BesselI(1, k*sigma)*b2*s1-B1*BesselK(1, k*sigma)*b2*s1+A1*BesselI(1, k*sigma)*b2+A1*BesselI(1, k*sigma)*sigma-B1*BesselK(1, k*sigma)*b2-B1*BesselK(1, k*sigma)*sigma)/sigma

> 

fn4 := T1*alpha1*(alpha1^2+T2)*(BesselK(1, alpha1)*b2*s1+BesselK(0, alpha1)*alpha1*b2+BesselK(1, alpha1)*b2-BesselK(1, alpha1))*c1-T1*alpha1*(alpha1^2+T2)*(BesselI(1, alpha1)*b2*s1-BesselI(0, alpha1)*alpha1*b2+BesselI(1, alpha1)*b2-BesselI(1, alpha1))*c2+T1*alpha2*(alpha2^2+T2)*(BesselK(1, alpha2)*b2*s1+BesselK(0, alpha2)*alpha2*b2+BesselK(1, alpha2)*b2-BesselK(1, alpha2))*c3-T1*alpha2*(alpha2^2+T2)*(BesselI(1, alpha2)*b2*s1-BesselI(0, alpha2)*alpha2*b2+BesselI(1, alpha2)*b2-BesselI(1, alpha2))*c4-T1*EZ*k*(K*k^2+K*T2+k1)*(A1*BesselI(1, k)*b2*s1-A1*BesselI(0, k)*b2*k-B1*BesselK(1, k)*b2*s1-B1*BesselK(0, k)*b2*k+A1*BesselI(1, k)*b2-B1*BesselK(1, k)*b2-A1*BesselI(1, k)+B1*BesselK(1, k));

T1*alpha1*(alpha1^2+T2)*(BesselK(1, alpha1)*b2*s1+BesselK(0, alpha1)*alpha1*b2+BesselK(1, alpha1)*b2-BesselK(1, alpha1))*c1-T1*alpha1*(alpha1^2+T2)*(BesselI(1, alpha1)*b2*s1-BesselI(0, alpha1)*alpha1*b2+BesselI(1, alpha1)*b2-BesselI(1, alpha1))*c2+T1*alpha2*(alpha2^2+T2)*(BesselK(1, alpha2)*b2*s1+BesselK(0, alpha2)*alpha2*b2+BesselK(1, alpha2)*b2-BesselK(1, alpha2))*c3-T1*alpha2*(alpha2^2+T2)*(BesselI(1, alpha2)*b2*s1-BesselI(0, alpha2)*alpha2*b2+BesselI(1, alpha2)*b2-BesselI(1, alpha2))*c4-T1*EZ*k*(K*k^2+K*T2+k1)*(A1*BesselI(1, k)*b2*s1-A1*BesselI(0, k)*b2*k-B1*BesselK(1, k)*b2*s1-B1*BesselK(0, k)*b2*k+A1*BesselI(1, k)*b2-B1*BesselK(1, k)*b2-A1*BesselI(1, k)+B1*BesselK(1, k))

Eqns are linear in c1,c2,c3,c4 - convert to matrix form

> 

A,b:=GenerateMatrix([fn1,fn2,fn3,fn4],[c1,c2,c3,c4]):

Solve the system with simpler variables

> 

AA:=Matrix(4,4,symbol=x):bb:=Vector(4,symbol=y):
> 

CC:=LinearSolve(AA,bb):

Substitute the more complicated expressions into the solution

> 

C:=eval(CC,{entries(AA=~A,'nolist'),entries(bb=~b,'nolist')}):
> 

soln:={c1=C[1],c2=C[2],c3=C[3],c4=C[4]}:

Check the solution

> 

simplify(eval([fn1,fn2,fn3,fn4],soln));

[0, 0, 0, 0]

NULL

Download 4eqns.mw

live connection...

dharr 8205
November 29 2023
2 1

I assume your tools->options->network has "enable maple cloud connection" checked? (becomes active next time you open Maple)

First 23 24 25 26 27 28 29 Last Page 25 of 81