dharr

Dr. David Harrington

8482 Reputation

22 Badges

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

Social Networks and Content at Maplesoft.com

Maple Application Center
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

Since text files are sequential and not random access, I think you have to live with this inefficiency, and just read and discard the first lines, say

filnam:=cat(currentdir(),"/graph8c.txt");
to 199 while readline(filnam) <> 0 do end do:
for i to 101 while (gr[i]:=readline(filnam))<>0 do end do:

and then ConvertGraph to get them in Maple Graph form. (The end of file conditions should probably be better handled here.)

You could read the whole file in at once with readdata or FileTools:-Text:-ReadFile if memory isn't a limitation.

I'd use map here:

map(x->`if`(x<0.2,0,x),a);

 

The index=real[2] is the problem. If this is changed to a regular index (index=1 in this case), then Minpoly works as expected.

restart;

alias(`~`[`=`](alpha__ || (1 .. 3), ` $`, RootOf(11*_Z^9+17*_Z^8-64*_Z^7-280*_Z^6+142*_Z^5+370*_Z^4+376*_Z^3-96*_Z^2+47*_Z-11, .2246 .. .2266), RootOf(11*_Z^9+17*_Z^8-64*_Z^7-280*_Z^6+142*_Z^5+370*_Z^4+376*_Z^3-96*_Z^2+47*_Z-11, 1.671 .. 1.68), RootOf(11*_Z^9+17*_Z^8-64*_Z^7-280*_Z^6+142*_Z^5+370*_Z^4+376*_Z^3-96*_Z^2+47*_Z-11, 2.648 .. 2.657)))

alpha__1, alpha__2, alpha__3

(1)

(Not all of these are independent; we could have chosen just one.)

evala(Algfield({`&alpha;__1`, `&alpha;__2`, `&alpha;__3`}))

({PDETools:-Solve})({`~`[`>=`](a, b, ` $`, 0), a^5*b+4*a^4*b^2+4*a^3*b^3-7*a^4*b-6*a^2*b^3-7*a*b^4+b^5-6*a^3*b+12*a^2*b^2+4*b^4+4*a^3-6*a*b^2+4*b^3+4*a^2-7*a*b+a = 0, a <> b})

{{a = alpha__1, b = RootOf(1216*_Z^4+(264*alpha__1^8+408*alpha__1^7-1580*alpha__1^6-6832*alpha__1^5+3508*alpha__1^4+9944*alpha__1^3+9948*alpha__1^2-10752*alpha__1+5204)*_Z^3+(891*alpha__1^8+1652*alpha__1^7-4748*alpha__1^6-24076*alpha__1^5+5354*alpha__1^4+35356*alpha__1^3+29668*alpha__1^2-196*alpha__1+3971)*_Z^2+(506*alpha__1^8+980*alpha__1^7-2264*alpha__1^6-12420*alpha__1^5+3676*alpha__1^4+11596*alpha__1^3+33800*alpha__1^2-7772*alpha__1+1210)*_Z-473*alpha__1^8-720*alpha__1^7+2560*alpha__1^6+10960*alpha__1^5-8034*alpha__1^4-13840*alpha__1^3-9304*alpha__1^2+1104*alpha__1-1133, index = real[2])}, {a = alpha__2, b = RootOf(1216*_Z^4+(264*alpha__2^8+408*alpha__2^7-1580*alpha__2^6-6832*alpha__2^5+3508*alpha__2^4+9944*alpha__2^3+9948*alpha__2^2-10752*alpha__2+5204)*_Z^3+(891*alpha__2^8+1652*alpha__2^7-4748*alpha__2^6-24076*alpha__2^5+5354*alpha__2^4+35356*alpha__2^3+29668*alpha__2^2-196*alpha__2+3971)*_Z^2+(506*alpha__2^8+980*alpha__2^7-2264*alpha__2^6-12420*alpha__2^5+3676*alpha__2^4+11596*alpha__2^3+33800*alpha__2^2-7772*alpha__2+1210)*_Z-473*alpha__2^8-720*alpha__2^7+2560*alpha__2^6+10960*alpha__2^5-8034*alpha__2^4-13840*alpha__2^3-9304*alpha__2^2+1104*alpha__2-1133, index = real[2])}, {a = alpha__3, b = RootOf(1216*_Z^4+(264*alpha__3^8+408*alpha__3^7-1580*alpha__3^6-6832*alpha__3^5+3508*alpha__3^4+9944*alpha__3^3+9948*alpha__3^2-10752*alpha__3+5204)*_Z^3+(891*alpha__3^8+1652*alpha__3^7-4748*alpha__3^6-24076*alpha__3^5+5354*alpha__3^4+35356*alpha__3^3+29668*alpha__3^2-196*alpha__3+3971)*_Z^2+(506*alpha__3^8+980*alpha__3^7-2264*alpha__3^6-12420*alpha__3^5+3676*alpha__3^4+11596*alpha__3^3+33800*alpha__3^2-7772*alpha__3+1210)*_Z-473*alpha__3^8-720*alpha__3^7+2560*alpha__3^6+10960*alpha__3^5-8034*alpha__3^4-13840*alpha__3^3-9304*alpha__3^2+1104*alpha__3-1133, index = real[2])}}

(2)

bSol := `~`[subs](%, b)

evalf[2*Digits](`~`[eval](11*_X^9-47*_X^8+96*_X^7-376*_X^6-370*_X^5-142*_X^4+280*_X^3+64*_X^2-17*_X-11, `~`[`=`](_X, bSol)))

{RootOf(1216*_Z^4+(264*alpha__1^8+408*alpha__1^7-1580*alpha__1^6-6832*alpha__1^5+3508*alpha__1^4+9944*alpha__1^3+9948*alpha__1^2-10752*alpha__1+5204)*_Z^3+(891*alpha__1^8+1652*alpha__1^7-4748*alpha__1^6-24076*alpha__1^5+5354*alpha__1^4+35356*alpha__1^3+29668*alpha__1^2-196*alpha__1+3971)*_Z^2+(506*alpha__1^8+980*alpha__1^7-2264*alpha__1^6-12420*alpha__1^5+3676*alpha__1^4+11596*alpha__1^3+33800*alpha__1^2-7772*alpha__1+1210)*_Z-473*alpha__1^8-720*alpha__1^7+2560*alpha__1^6+10960*alpha__1^5-8034*alpha__1^4-13840*alpha__1^3-9304*alpha__1^2+1104*alpha__1-1133, index = real[2]), RootOf(1216*_Z^4+(264*alpha__2^8+408*alpha__2^7-1580*alpha__2^6-6832*alpha__2^5+3508*alpha__2^4+9944*alpha__2^3+9948*alpha__2^2-10752*alpha__2+5204)*_Z^3+(891*alpha__2^8+1652*alpha__2^7-4748*alpha__2^6-24076*alpha__2^5+5354*alpha__2^4+35356*alpha__2^3+29668*alpha__2^2-196*alpha__2+3971)*_Z^2+(506*alpha__2^8+980*alpha__2^7-2264*alpha__2^6-12420*alpha__2^5+3676*alpha__2^4+11596*alpha__2^3+33800*alpha__2^2-7772*alpha__2+1210)*_Z-473*alpha__2^8-720*alpha__2^7+2560*alpha__2^6+10960*alpha__2^5-8034*alpha__2^4-13840*alpha__2^3-9304*alpha__2^2+1104*alpha__2-1133, index = real[2]), RootOf(1216*_Z^4+(264*alpha__3^8+408*alpha__3^7-1580*alpha__3^6-6832*alpha__3^5+3508*alpha__3^4+9944*alpha__3^3+9948*alpha__3^2-10752*alpha__3+5204)*_Z^3+(891*alpha__3^8+1652*alpha__3^7-4748*alpha__3^6-24076*alpha__3^5+5354*alpha__3^4+35356*alpha__3^3+29668*alpha__3^2-196*alpha__3+3971)*_Z^2+(506*alpha__3^8+980*alpha__3^7-2264*alpha__3^6-12420*alpha__3^5+3676*alpha__3^4+11596*alpha__3^3+33800*alpha__3^2-7772*alpha__3+1210)*_Z-473*alpha__3^8-720*alpha__3^7+2560*alpha__3^6+10960*alpha__3^5-8034*alpha__3^4-13840*alpha__3^3-9304*alpha__3^2+1104*alpha__3-1133, index = real[2])}

 

{-0.7765721e-11, -0.40e-16, -0.2e-17}

(3)

evalf(bSol)

{.3768363389, .5974151794, 4.441922439}

(4)

Again, these may not be independent - this time Algfield fails

PA := evala(Algfield(bSol)); PA[4]

false

(5)

Could this be because of the real[2] indexing? Find the simple index corresponding to it. In each case it is index=1

evalf([seq(subs(real[2] = i, bSol[1]), i = 1 .. 4)]); evalf([seq(subs(real[2] = i, bSol[2]), i = 1 .. 4)]); evalf([seq(subs(real[2] = i, bSol[3]), i = 1 .. 4)])

[HFloat(0.3768363389232699), -HFloat(1.1781382393831348)+HFloat(1.6342731836523705)*I, HFloat(-0.8215019563741058), -HFloat(1.1781382393831348)-HFloat(1.6342731836523705)*I]

 

[HFloat(4.44192244093814), HFloat(0.7521012471815752)+HFloat(0.18644586398558147)*I, HFloat(-2.670902846485502), HFloat(0.7521012471815752)-HFloat(0.18644586398558147)*I]

 

[HFloat(0.5974151764626043), HFloat(8.500619958633818)+HFloat(3.0727707478333337)*I, HFloat(-3.6203669440591857), HFloat(8.500619958633818)-HFloat(3.0727707478333337)*I]

(6)

bSol2 := subs(real[2] = 1, bSol)

{RootOf(1216*_Z^4+(264*alpha__1^8+408*alpha__1^7-1580*alpha__1^6-6832*alpha__1^5+3508*alpha__1^4+9944*alpha__1^3+9948*alpha__1^2-10752*alpha__1+5204)*_Z^3+(891*alpha__1^8+1652*alpha__1^7-4748*alpha__1^6-24076*alpha__1^5+5354*alpha__1^4+35356*alpha__1^3+29668*alpha__1^2-196*alpha__1+3971)*_Z^2+(506*alpha__1^8+980*alpha__1^7-2264*alpha__1^6-12420*alpha__1^5+3676*alpha__1^4+11596*alpha__1^3+33800*alpha__1^2-7772*alpha__1+1210)*_Z-473*alpha__1^8-720*alpha__1^7+2560*alpha__1^6+10960*alpha__1^5-8034*alpha__1^4-13840*alpha__1^3-9304*alpha__1^2+1104*alpha__1-1133, index = 1), RootOf(1216*_Z^4+(264*alpha__2^8+408*alpha__2^7-1580*alpha__2^6-6832*alpha__2^5+3508*alpha__2^4+9944*alpha__2^3+9948*alpha__2^2-10752*alpha__2+5204)*_Z^3+(891*alpha__2^8+1652*alpha__2^7-4748*alpha__2^6-24076*alpha__2^5+5354*alpha__2^4+35356*alpha__2^3+29668*alpha__2^2-196*alpha__2+3971)*_Z^2+(506*alpha__2^8+980*alpha__2^7-2264*alpha__2^6-12420*alpha__2^5+3676*alpha__2^4+11596*alpha__2^3+33800*alpha__2^2-7772*alpha__2+1210)*_Z-473*alpha__2^8-720*alpha__2^7+2560*alpha__2^6+10960*alpha__2^5-8034*alpha__2^4-13840*alpha__2^3-9304*alpha__2^2+1104*alpha__2-1133, index = 1), RootOf(1216*_Z^4+(264*alpha__3^8+408*alpha__3^7-1580*alpha__3^6-6832*alpha__3^5+3508*alpha__3^4+9944*alpha__3^3+9948*alpha__3^2-10752*alpha__3+5204)*_Z^3+(891*alpha__3^8+1652*alpha__3^7-4748*alpha__3^6-24076*alpha__3^5+5354*alpha__3^4+35356*alpha__3^3+29668*alpha__3^2-196*alpha__3+3971)*_Z^2+(506*alpha__3^8+980*alpha__3^7-2264*alpha__3^6-12420*alpha__3^5+3676*alpha__3^4+11596*alpha__3^3+33800*alpha__3^2-7772*alpha__3+1210)*_Z-473*alpha__3^8-720*alpha__3^7+2560*alpha__3^6+10960*alpha__3^5-8034*alpha__3^4-13840*alpha__3^3-9304*alpha__3^2+1104*alpha__3-1133, index = 1)}

(7)

This time Algfield succeeds, and we only need the single algebraic number alpha__1

PA := evala(Algfield(bSol2)); PA[4]; K := PA[3]

true

 

{RootOf(11*_Z^9+17*_Z^8-64*_Z^7-280*_Z^6+142*_Z^5+370*_Z^4+376*_Z^3-96*_Z^2+47*_Z-11, index = 1)}

(8)

`~`[`@`(evala, Minpoly)](bSol2, _X)

{_X^9-(47/11)*_X^8+(96/11)*_X^7-(376/11)*_X^6-(370/11)*_X^5-(142/11)*_X^4+(280/11)*_X^3+(64/11)*_X^2-(17/11)*_X-1}

(9)

`~`[PolynomialTools[MinimalPolynomial]](bSol2, _X)

{11*_X^9-47*_X^8+96*_X^7-376*_X^6-370*_X^5-142*_X^4+280*_X^3+64*_X^2-17*_X-11}

(10)

NULL

Download minpoly.mw

It is a sum over the different roots. Perhaps there is an easier way, but here is one way to do it.

restart

S := proc (t) options operator, arrow; tanh(ln(1+t^2)) end proc

proc (t) options operator, arrow; tanh(ln(1+t^2)) end proc

q := diff(S(t), `$`(t, n))

pochhammer(1-n, n)+Sum(-(1/4)*_alpha^3*pochhammer(-n, n)*(t-_alpha)^(-1-n), _alpha = RootOf(_Z^4+2*_Z^2+2))

`assuming`([simplify(allvalues(q))], [n::posint])

(1/4)*GAMMA(n+1)*(-1)^n*((1+I)*(-1-I)^(1/2)*(t-(-1-I)^(1/2))^(-1-n)+(1-I)*(-1+I)^(1/2)*(t-(-1+I)^(1/2))^(-1-n)+(-1-I)*(-1-I)^(1/2)*(t+(-1-I)^(1/2))^(-1-n)+(-1+I)*(t+(-1+I)^(1/2))^(-1-n)*(-1+I)^(1/2))

``

Download simplify.mw

Your code contains the line

S2 := t -> piecewise(xn(t) <= xa and 0 < zn(t), -S1(t)*exp(mu1*alpha2(t)), 0)

which suggests that the drop-off might be when xn(t) or zn(t) becomes zero. Plotting xn(t) or equivalently, Xn, shows that Xn goes through zero at that point. So then fsolve(Xn, 0.2..0.9) returns the drop-off point of 0.5006805969.

min_problem.mw

Use Pi and not pi.

restart; with(LinearAlgebra)

w := c[i]*(1-cos(2*Pi*x/a))*(1-cos(2*Pi*y/b))

c[i]*(1-cos(2*Pi*x/a))*(1-cos(2*Pi*y/b))

al_eq1 := (1/2)*(int(int(D__11*(diff(w, x, x))^2+2*D__12*(diff(w, x, x))*(diff(w, y, y))+4*D__66*(diff(w, x, y))^2+D__22*(diff(w, y, y))^2-2*q__0*w, x = 0 .. a), y = 0 .. b))

c[i]*(6*D__11*Pi^4*b^5*c[i]+4*D__12*Pi^4*a^2*b^3*c[i]+6*D__22*Pi^4*a^4*b*c[i]+8*D__66*Pi^4*a^2*b^3*c[i]-a^4*b^5*q__0)/(a^3*b^4)

NULL

Download Pi.mw

A bit less efficient, perhaps.

[Edit: my interpretation of the OPs question was that different columns were to be scaled differently, but if rows were intended, then premultiplication by the diagonal matrix works.]

A:=Matrix(3,3,[1,6,9,7,4,3,2,8,9]);
B:=<a,b,c>;
A^~(-1).LinearAlgebra:-DiagonalMatrix(B);

Matrix(3, 3, {(1, 1) = 1, (1, 2) = 6, (1, 3) = 9, (2, 1) = 7, (2, 2) = 4, (2, 3) = 3, (3, 1) = 2, (3, 2) = 8, (3, 3) = 9})

Vector(3, {(1) = a, (2) = b, (3) = c})

Matrix(%id = 36893490033002534060)

``

Download div.mw

If I understand correctly what you want, a (Tabulated) DataFrame is one way to display this.

stiffnessmatrixcal.mw

Here's one way to do it, using custom tickmarks.

restart;

Frequency as a function of lambda

freq:=lam->2.9979/lam;

proc (lam) options operator, arrow; 2.9979/lam end proc

fmax:=0.3;

.3

lambda values for ticks

lams:=[10,15,20,40];
freqs:=map(freq,lams);

[10, 15, 20, 40]

[.2997900000, .1998600000, .1498950000, 0.7494750000e-1]

p1:=plot(x*(1-x),x=0..1,view=[default,0..fmax],labels=[x,f],axes=boxed):
p2:=plot(x*(1-x),x=0..1,view=[default,0..fmax],labels=[x,lambda],axes=boxed,axis[2]=[tickmarks=(freqs=~lams)]):
plots:-dualaxisplot(p1,p2);

NULL

Download dualaxes.mw

Suggest you upload your worksheet next time (green up arrow in the editor).

in
for n to N do

subs([seq(x(i) = x[c][i], i = 0 .. n),

the i in x[c][i] must be x[c][i+1] since list indices start at 1. Same for further on in the line. After fixing this there is another error arising because you have used "x" as a simple variable in the differentiation and other places as well as with x[c][i] - you need different names for these. 

There are various technical issues here that I'm going to ignore. What answer were you expecting? What do you want to happen for negative x? Is this a signal that is zero for x<0. If so, then you can just use laplace with s=I*omega

eval(laplace(f(x),x,s),s=I*omega);

gives
 1/2*Pi^(1/2)/(I*omega)^(3/2)

Well, this finds any of the polynomials. But what type of answer exactly is expected here?

restart

eq := (-x^2+x)^(4*n) = (x^2+1)*p+(-4)^n

(-x^2+x)^(4*n) = (x^2+1)*p+(-4)^n

p := solve(eq, p)

-((-4)^n-(-x^2+x)^(4*n))/(x^2+1)

p7 := normal(eval(p, n = 7), expanded)

x^54-28*x^53+377*x^52-3248*x^51+20098*x^50-95032*x^49+356642*x^48-1089008*x^47+2751463*x^46-5817892*x^45+10371647*x^44-15656288*x^43+20050108*x^42-21785872*x^41+20066492*x^40-15656288*x^39+10355263*x^38-5817892*x^37+2767847*x^36-1089008*x^35+340258*x^34-95032*x^33+36482*x^32-3248*x^31-16007*x^30-28*x^29+16385*x^28-16384*x^26+16384*x^24-16384*x^22+16384*x^20-16384*x^18+16384*x^16-16384*x^14+16384*x^12-16384*x^10+16384*x^8-16384*x^6+16384*x^4-16384*x^2+16384

NULL

Download polynomial.mw

use sol:-value(r=0..10,t=0..20) to get the values out. [Edit - not correct - see correct answer below] (It is helpful if you load your actual worksheet using the big green up-arrow in the Mapleprimes editor.)

The error is because C1 is zero. My recollection is that the calculation of the t-value involves dividing by the standard error, which is zero. Then the absolute value of infinity leads to the complex infinity.

You should scale your dependent data down to values within a few orders of magnitude of 1, say by dividing it by 10^7. Then your parameters will be 10^7 times higher and more reasonable. You are likely running into numerical accuracy issues because your problem is not better scaled. Actually you could use stronger scaling for the last column.

Here is a scaled version that works:

toycode3.mw

with(GraphTheory);
G := Multigraph(2, [{1, 2}, {1, 2}]);
DrawGraph(G);

The 2 by the edge in DrawGraph shows there are two edges between vertices 1 and 2.

First 32 33 34 35 36 37 38 Last Page 34 of 83