mmcdara

6730 Reputation

18 Badges

8 years, 181 days

MaplePrimes Activity


These are replies submitted by mmcdara

@michaelvio 

You've put my mind at rest because I really felt we were going round in circles and I didn't understand what you were saying.

@felixiao 

In my previous code I use notations like 'k0' to prevent the replacement of the name k0 by its numeric value and thus ease the comprison of the two matrices.

If you want to eval the matrices remove the single quotes:
qqs_mmcdara_0.mw
qqs_mmcdara_1.mw
qqs_mmcdara_2.mw

@michaelvio 

J := x0 -> evalf(Int(T, x = x0 .. x0+1.542976947e-13));

J(1/528);
J(1/435);

# and so on

If it is not what you want find someone else, I'm done

@michaelvio 

Firstly the integrand is far too complex to expect that int may return a value: use evalf@Int instead (see ?evalf/Int for more details).

But even with that the values of the integrand are far beyond the limitations of Maple to get any accute estimation of the integral:

restart

a := -1.44670357887361*10^(-7);

-0.1446703579e-6

 

-0.1049267156e-8

 

0.1890440485e-11

 

-0.6233924848e-15

 

0.762014687e-2*t-0.1446703579e-6*t^2-0.1049267156e-8*t^3+0.1890440485e-11*t^4-0.6233924848e-15*t^5

 

0.762014687e-2-0.2893407158e-6*t-0.3147801468e-8*t^2+0.7561761940e-11*t^3-0.3116962424e-14*t^4

 

3.259000000

 

3.950000001

 

3.2590001

(1)

T := subs(t = 1/x, 5.012764943*10^(-24)*Ea/(exp(Ea/(4.100527530*10^(-21)))-1))/x^2;

0.5012764943e-23*(0.762014687e-2/x-0.1446703579e-6/x^2-0.1049267156e-8/x^3+0.1890440485e-11/x^4-0.6233924848e-15/x^5)/((exp(0.1858333303e19/x-0.3528091369e14/x^2-0.2558858947e12/x^3+461023727.0/x^4-152027.3867/x^5)-1)*x^2)

(2)

J := Int(T, x = 0.2298850575e-2 .. 0.2298850575e-2+1.542976947*10^(-13));

Int(0.5012764943e-23*(0.762014687e-2/x-0.1446703579e-6/x^2-0.1049267156e-8/x^3+0.1890440485e-11/x^4-0.6233924848e-15/x^5)/((exp(0.1858333303e19/x-0.3528091369e14/x^2-0.2558858947e12/x^3+461023727.0/x^4-152027.3867/x^5)-1)*x^2), x = 0.2298850575e-2 .. 0.2298850575e-2)

(3)

evalf(J);

# as the upper bound is extremely closed to the upper bound, it is likely that
# the formula below gives an acceptable approximation of J...unlesss T is
# rapidly cha,ging in the range 0.2298850575e-2 .. 0.2298850575e-2+1.542976947*10^(-13)
# (which is not the case as x=0.2298850575e-2 is not a root of denom(T))=


'J' = 'eval(T, x = 0.2298850575e-2)' * 1.542976947*10^(-13)

0.

 

J = 0.1542976947e-12*(eval(T, x = 0.2298850575e-2))

(4)

# Note that x=0.002298850575 is not a root of (numer(T).
# Look at this a priori strange result

eval(T, x=0.002298850575);

0.

(5)

# Why this result?

TND := simplify~(normal(T));

NUM := eval(numer(TND), x=0.002298850575);
DEN := eval(denom(TND), x=0.002298850575);  

# Explanation

UnderExponential := op([1, 1, 1], denom(TND));

UnderExponential := eval(UnderExponential, x=0.002298850575);

(0.3819800509e-25*x^4-0.7251984984e-30*x^3-0.5259729615e-32*x^2+0.9476333790e-35*x-0.3124919994e-38)/((exp((0.1858333303e19*x^4-0.3528091369e14*x^3-0.2558858947e12*x^2+461023727.*x-152027.3867)/x^5)-1.)*x^7)

 

0.1048854911e-35

 

Float(infinity)

 

(0.1858333303e19*x^4-0.3528091369e14*x^3-0.2558858947e12*x^2+461023727.*x-152027.3867)/x^5

 

0.7947757879e21

(6)

Maple_floats(MAX_FLOAT);

0.1e9223372036854775807

(7)

# So exp(UnderExponential) is a number with more than 10^20 digits
# As this numer exceeds

Maple_floats(MAX_FLOAT);

0.1e9223372036854775807

(8)

# Here is the minimum value (X) of x wuch that exp(UnderExponential) is notundefined:

UnderExponential := op([1, 1, 1], denom(TND)):

X := 'X':
X := fsolve(eval(UnderExponential, x=X) / log(10) = op(2, Maple_floats(MAX_FLOAT)));

exp(eval(UnderExponential, x=X));

print():
exp(eval(UnderExponential, x=X - 1e-10));

0.8748151149e-1

 

0.2996941173e9223372028516350198

 

 

Float(infinity)

(9)
 

 

Download AphS_mmcdara.mw

@GFY 

fsolve can provide several solutions see ?fsolve/details (option maxsols).

@Kitonum 

Good you answered, I remove my comment

@michaelvio 

I believe the idea was not that bad but I should have been more cautious and think twice before acting.

The minimization does not return an acceptable value for A.
The residue E1-F(Ea) is equal to E1 because F(Ea) is almost equal to 0 at the optimum A.

I belive this problem comes from the fact E1 is doubly defined: either as a the derivative of Ea,or either by the relation E1=F(Ea)
This leads to an ODE Ea must verify... but the expession 

Ea := 0.762014687e-2*t+a*t^2+b*t^3+c*t^4+d*t^5

is not a solution of this ode.

I would say that something is not correct in your problem, see here

restart

0.762014687e-2*t+a*t^2+b*t^3+c*t^4+d*t^5

 

0.762014687e-2+2*a*t+3*b*t^2+4*c*t^3+5*d*t^4

 

3.314763888+189225*a+82312875*b+35806100625*c+15575653771875*d

 

4.023437547+278784*a+147197952*b+77720518656*c+41036433850368*d

 

21.79362005+8179600*a+23393656000*b+66905856160000*c+191350748617600000*d

(1)

# What you do later is
#
# solve(
#   {E1 = (5.012764943*10^(-24))*Ea/(exp(Ea/(4.100527530*10^(-21)))-1), E2 = 3.259, E3 = 3.95, E4 = 3.259}
#   , {a, b, c, d}
# )
#
# Let's focus on the first relation.
# It writes more simply E1 = F(Ea) where F is the function defined by
# (p = 5.012764943*10^(-24) and q = 4.100527530*10^(-21)

F := Z -> (p*Z*exp(-Z/q)-1)

proc (Z) options operator, arrow; p*Z*exp(-Z/q)-1 end proc

(2)

# Let's put Z = Ea(t).
#
# This means that E1(t) has two definitions:

`<|>`(`either :`, 'E1'(t) = Diff('Ea'(t), t));

`<|>`(`or :`, 'E1'(t) = 'F'('Ea'(t)));

`<|>`(`and thus :`, Diff('Ea'(t), t) = 'F'('Ea'(t)));

`<|>`(`or explicitely`, Diff('Ea'(t), t) = F('Ea'(t)));

Matrix([[`either :`, E1(t) = Diff(Ea(t), t)]])

 

Matrix([[`or :`, E1(t) = F(Ea(t))]])

 

Matrix([[`and thus :`, Diff(Ea(t), t) = F(Ea(t))]])

 

Matrix([[`or explicitely`, Diff(Ea(t), t) = p*Ea(t)*exp(-Ea(t)/q)-1]])

(3)

# Solution of the previous ode:

value((3)[1, 2]);

dsolve(%);

diff(Ea(t), t) = p*Ea(t)*exp(-Ea(t)/q)-1

 

t-Intat(1/(p*_a*exp(-_a/q)-1), _a = Ea(t))+_C1 = 0

(4)

# Replacing Ea(t) by omega*t+a*t^2+b*t^3+c*t^4+d*t^5 and setting t=0 shows
# omega should beequal to -1 to expect that the relation could be satistied
# for all t >= 0.

value((3)[1, 2]);

EQ := eval(%, Ea(t) = omega*t+a*t^2+b*t^3+c*t^4+d*t^5);

eval(EQ, t=0);

diff(Ea(t), t) = p*Ea(t)*exp(-Ea(t)/q)-1

 

5*d*t^4+4*c*t^3+3*b*t^2+2*a*t+omega = p*(d*t^5+c*t^4+b*t^3+a*t^2+omega*t)*exp(-(d*t^5+c*t^4+b*t^3+a*t^2+omega*t)/q)-1

 

omega = -1

(5)
 

 

Download Aph1_Pb_0.mw


Firstly couple is never defined in your worksheet; I guess you meant totalsecular2 instead?

Based on the this assumption I propose you this question929_mmcdara.mw.

It's not clear whether you want to use solve or fsolve?
(note that the numeric values in totalsecular2 being floats, solve will return solutions with floats, but will behave as usual, meaning it will be able to deliver several solutions).

You will see there are no real solutions, so fsolve requires the option complex to find one


Here is an add-on.mw where I replaced sin(varsigma) by s and cos(varsigma) by sqrt(1-s^2).Ypu wil see that the number of solutions solve returns can be much lower while they still remain complex.

@grompvevo 

Shouldn't you speak of weak from and strong form of conservation laws instead of global and local forms?
Weak and Strong forms

@felixiao 

Sorry I missed your file in your initial question

restart
for r to 3 do 
  T[r] := exp(r+1); 
  P[r] := exp(r-1) 
end do:

mul(T[r]*P[r], r=1..3)
                         2                     
                 (exp(2))  exp(3) exp(1) exp(4)

or

restart
TP := 1:

for r to 3 do 
  T[r] := exp(r+1); 
  TP   := TP*T[r]:
  P[r] := exp(r-1):
  TP   := TP*P[r]:
end do:

TP
                         2                     
                 (exp(2))  exp(3) exp(1) exp(4)

 

@acer 

The equation the OP presents is the one dimensional (N=1) wave equation. This equation derive from a variational principle of energy conservation.


A good way to test a numerical scheme is to verify that E(t) is indeed a constant when u(x, t) is replaced by a numerical approximation uN(x, t)
This is how I would interpret "I want to test how well Maple's numeric solver respects a variety of conservation laws"

As a rule stable numerical schemes dissipate (mass, impulse, energy) and E(t) is likely do decrease over time. The smaller this decrease the better the solving method.
From your file:

restart;

pde := diff(v(x, t), t, t) = diff(v(x, t), x, x):

f := xi -> exp(-xi^2):

a := -10: b := 10: dx := 1/50: t_final := 10:

pds := pdsolve(pde, {v(a, t) = 0, v(b, t) = 0, v(x, 0) = f(x + 5),
                     D[2](v)(x, 0) = -eval(diff(f(x), x), x = x + 5)},
               numeric, range = a .. b, time = t, spacestep = dx):

sol_proc := rhs(pds:-value(output = listprocedure)[3]):

sol := (x, t) -> piecewise(a < x and x < b, sol_proc(x, t), 0):

# plot(D[2](sol)(x,0), x=a..b, size=[500,200]);

E := t ->
     evalf(Int((D[2](sol)(x,t))^2 + (D[1](sol)(x,t))^2, x=a..b, digits=10, epsilon=1e-5, method=_d01ajc))    ;

proc (t) options operator, arrow; evalf(Int((D[2](sol))(x, t)^2+(D[1](sol))(x, t)^2, x = a .. b, digits = 10, epsilon = 0.1e-4, method = _d01ajc)) end proc

(1)

E(0)

2.507375287

(2)

for s from 1 to 5 do
  [s, E(s)]
end do;

[1, 2.505877784]

 

[2, 2.505877752]

 

[3, 2.505877729]

 

[4, 2.505877918]

 

[5, 2.505873355]

(3)
 

 

Download E.mw



With Maple 2015:

M21 := expand(M21);
                               0

M_final := MatrixInverse(Matrix(2, 2, {(1, 1) = 1/M11, (1, 2) = -M12/M11, (2, 1) = 0, (2, 2) = -M22/M21}));
%;
Error, (in MatrixInverse) numeric exception: division by zero


 

@JAMET 

Pourquoi n'avez-vous pas dit plus tôt que les codes qui accompagnent vos questions sont générés par une IA (à l'évidence pas si intelligente que ça)?

Cela aurait peut-être permis de comprendre pourquoi ils sont un tel foutoir et pourquoi vos questions successives laissent l'impression que vous n'avez rien appris, au fil des ans, des précédentes réponses que vous avez reçues.
(Vous êtes vous aperçu que le texte que vous avez fourni dans votre question contient deux séquences débitant chacune par un restart ? N'auriez-vous pas pu le nettoyer un peu avant de le livrer afin qu'il soit plus présentable ?),

En outre pourquoi ne délivrez-vous jamais de feuille de travail Maple ?
Si vous n'avez pas de licence personnelle (ce que je peux comprendre car elle est payante) à quoi vous servent alors les codes contenus dans les réponses que vous recevez?
Un petit éclaircissement serait le bienvenu.

Si vous êtes passionné de géométrie, pourquoi ne pas utiliser Geogebra qui est gratuit et bénéficie d'un excellent support (en français)... et qui en outre propose es applications permettant de tracer les hyperbole d'Apollonius d'une conique (dont je n'ai trouvé aucune trace par ailleurs) ?

Pour conclure cette réponse, et je vous prie une fois de plus d'excuser mes propos dans ma réponse précédente, avoir 85 ans et être defficient visuel ne devrait pas être un prétexte pour ne pas respecter les quelques règles de courtoisie en vigueur sur ce site.


(improved DEEPL translation at the end of the comment)

 restart;

with(geometry):
_EnvHorizontalName := x: _EnvVerticalName := y:

F := proc(__A, __B, __C, ratio)

  local A, B, C, T, G, XY, AB, BC, CA, inC, exC, lAB, lBC, lCA;
  local H1AB, H1CA, H2CA, H2BC, H3AB, H3BC, H:
  local f, Apo, tx, TX, i, P, GP, ici;

  uses geometry, plots:
  
  point(A, op(__A)):
  point(B, op(__B)):
  point(C, op(__C)):

  triangle(T, [A, B, C]);
  centroid(G, T);
  XY := coordinates(G);

  bisector(AB, C, T);
  bisector(BC, A, T);
  bisector(CA, B, T);   
  
  incircle(inC, T);
  excircle(exC, T, [c1(o1),c2(o2),c3(o3)]);

  line(lAB, [A, B]):
  line(lBC, [B, C]):
  line(lCA, [C, A]):

  f   := (x, y) -> sqrt((x - __A[1])^2 + (y - __A[2])^2)/sqrt((x - __B[1])^2 + (y - __B[2])^2) - ratio;
  Apo := implicitplot(
           f(x, y)
           , x = -5 .. 5, y = -5 .. 5
           , grid = [100, 100]
           , style = line
           , color = red
           , thickness = 2
          );

  tx := NULL:

  TX := ["A", "B", "C"]:
  for i from 1 to 3 do
    P   := args[i]:
    GP  := P -~ XY:
    ici := P +~ GP*~0.1:
    tx := tx, textplot([ici[],  TX[i]], font = [times, bold, 16])
  end do:

intersection(H1AB, exC[1], lAB);
intersection(H1CA, exC[1], lCA);
intersection(H2CA, exC[2], lCA);
intersection(H2BC, exC[2], lBC);
intersection(H3AB, exC[3], lAB);
intersection(H3BC, exC[3], lBC);
H := evalf(coordinates~([H1AB, H1CA, H2CA, H2BC, H3AB, H3BC, A, B, C]));

  plots:-display(
    [
      Apo
      , tx
      , draw(T, color = black, thickness = 2)
      , draw([AB, BC, CA], style = line, color = blue, linestyle=3)
      , draw(inC, style = line, color = green, thickness = 2)
      , draw(exC, style = line, color = magenta, thickness = 2)
      , draw([lAB, lBC, lCA], color=black, linestyle=4)
    ]
    , axes = none
    , view=[((min-1)..(max+1))(map2(op, 1, H)), ((min-1)..(max+1))(map2(op, 2, H))]
    , size=[600, 600]
    , scaling = constrained
    , title = "Triangle with Bisectors, Apollonius Hyperbola, and Circles"
  );

end proc:

A := [0, 0];

B := [4, 2];

C := [2, 3];
F(A, B, C, 3)

 

[0, 0]

 

[4, 2]

 

[2, 3]

 

 

 

 

Download Proposition.mw


 

Why didn't you say earlier that the text codes accompanying  your questions are generated by an AI (obviously not that intelligent)?

That might have helped us to understand why they're such a mess, and why your repeated questions leave the impression that you've learned nothing, over the years, from the previous answers you've received.
(Did you notice that the text you provided in your last question contains two sequences, each beginning with a restart? Couldn't you have cleaned it up a bit before delivering and make it more presentable?),

And besides, why don't you ever deliver any Maple worksheets?
If you don't have a personal license (which I can understand, as it's not free), then what use are the codes uploaded in the answers you receive?
A little clarification would be welcome.

If you're passionate about geometry, why not use Geogebra, which is free and benefits from an excellent french support... and which also offers free applications for tracing Apollonius hyperbolas of a conic (of which I've found no mention elsewhere)?

To conclude this reply, and once again I apologize for my words in my previous comment, being 85 years old and visually defficient should not be an excuse for not respecting the few rules of courtesy that apply on this site.

Translated with DeepL.com (free version)

 

@Carl Love 

It works perfectly well and goes far beyond what I expected.

You told me once that subsindets was one of the most powerful algebraic operation in Maple. But I'm still not at ease with it and keep using makeshift alternatives based on things like map and select, or patmatch sometimes.
Since I read your last reply I've been practicing rewriting a few things using subsindets (whose help page is not very clear, in particular concerning the construction of fcntype), and it's very powerful indeed. 
Now all I have to do is give up my old habits.

Thanks for your time.

@nm 

Your reply suits me.
I can consider that I found a solution to what to a limitation of my old 2015 version and that it has been fixed in the most recent versions. 
I can't decently ask someone to explain me why something which didn't work 10 years ago works today: codes evolve and it's all for the better.

Thanks again

By the way: in Maple 2015 you can convert an expression containing intat/Intat into another one containing int/Int, but the reciprocal is not possible. Is that possible in Maple 2024?
For instance:

kernelopts(version)

`Maple 2015.2, APPLE UNIVERSAL OSX, Dec 20 2015, Build ID 1097895`

(1)

J := Intat(f(x), x=y)

Intat(f(x), x = y)

(2)

K := convert(J, Int);

Int(f(y), y)

(3)

L := convert(K, Intat)

Error, unrecognized conversion

 

# one might have expected something like

L = Intat(f(_Z), _Z=y)

L = Intat(f(_Z), _Z = y)

(4)
 

 

Download Intat.mw

2 3 4 5 6 7 8 Last Page 4 of 142