mmcdara

A question ofsize....

Answered: mmcdara 6149

April 23 2023

2 6

A remark: had you download your worksheet by using the green up-arrow in the menu bar, it would have been easier for anyone here to work on your problem. Think to it the next time you post a question :-)

As I said in the title the system is formally simple (a 4-by-for linear system of full rank) but the expressions of the coefficients of this system are that complex that solve seems incapable to solve it without a little help.

You will find here a simple idea: solve a formal 4-by-4 linear system and, AFTER, replace the coefficients of the solution you got by the true coefficients.

>

restart

>	interface(version)

(1)

>

u := c1*BesselK(0, alpha1*r) + c2*BesselI(0, alpha1*r) + c3*BesselK(0, alpha2*r) + c4*BesselI(0, alpha2*r) - A1*k^2*(k^2 - (-alpha^2*j*s + 4*c*s)/c)*BesselK(0, k*r)/((-alpha1^2 + k^2)*(-alpha2^2 + k^2)) - B1*k^2*(k^2 - (-alpha^2*j*s + 4*c*s)/c)*BesselI(0, k*r)/((-alpha1^2 + k^2)*(-alpha2^2 + k^2)) + ur*(-alpha^2*j*s + 4*c*s)/(c*(1 + c)):

>	w := (-1 - c)/(2s(4.c - alpha^2j))collect(diff(diff(rdiff(u, r), r)/r, r) + (alpha^2 - 4cs)diff(u, r)/(1 + c) - A1k^3BesselK(1, kr) + B1k^3BesselI(1, k*r), [c1, c2, c3, c4], factor):

>	om := (-1)/2*diff(u, r)/u:

>	fn1 := collect(simplify(subs(subs(r = 1, u))), [c1, c2, c3, c4], factor):

>	fn2 := collect(simplify(subs(subs(r = sigma, u))), [c1, c2, c3, c4], factor):

>	fn3 := collect(simplify(subs(r = 1, w)), [c1, c2, c3, c4], factor):

>	fn4 := collect(simplify(subs(r = sigma, w)), [c1, c2, c3, c4], factor):

>	DesiredSolution := solve({fn1 = 0, fn2 = 0, fn3 = 0, fn4 = 0}, {c1, c2, c3, c4})

Warning, computation interrupted

>	# Check the system Vector[row](4, i -> has(fn\|\|i, [c1, c2, c3, c4])); Matrix(4$2, (i, j) -> type(fn\|\|i, linear([c\|\|j])));

Matrix([[true, true, true, true], [true, true, true, true], [true, true, true, true], [true, true, true, true]])

(2)

>	# This because I'm going to chose two names that are not already used indets({fn1, fn2, fn3, fn4}, name)

(3)

>

# This is a 4-by-4 linear system wrt [c1, c2, c3, c4].
#
# Let MatrixCoeffs be the matrix of coefficients of c(j) in fn(i).
# Let RhsCoeffs bethe column vector of terms in fn(i) which do not
# contain c1, ..., c4.

MatrixCoeffs := Matrix(4$2, (i, j) -> coeff(fn||i, c||j)):
RhsCoeffs := Vector(4, i -> -eval(fn||i, [c1, c2, c3, c4]=~0)):

>	# Solve this symbolic system: SymbolicMatrix := Matrix(4$2, symbol=m): SymbolicRhs := Vector(4, symbol=r): SymbolicSol := LinearAlgebra:-LinearSolve(SymbolicMatrix, SymbolicRhs): SymbolicSol := simplify(SymbolicSol, size): # Example SymbolicSol[1]

(4)

>	# Here is the lengthof SymbolicSol: length(SymbolicSol)

(5)

>

# So try to imagine what would be the length of Desiredsolution... which is
# the reason why solve faces huge problems.
#
# But, knowing SymbolicSol, MatrixCoeffs and RhsCoeffs, it becomes extremely
# simple to get the expression of any c(j).

Substitutions := {entries(SymbolicMatrix =~ MatrixCoeffs, nolist), entries(SymbolicRhs =~ RhsCoeffs, nolist)}:

C_Solution := Vector(4, j -> eval(SymbolicSol[j], Substitutions)):

>	# For instance length(C_Solution[1]); C1 := C_Solution[1]:

(6)

>

Download LinearSolve.mw

Simplifying the expressions of any element of C_Solution is now a very complex task (if possible).

You have two different data, let's s...

Answered: mmcdara 6149

April 21 2023

0 0

You have two different data, let's say B and Y, even if Y never appears here, and you knox, or suspect,there is among them a functional relation Y = a*B+b*sinh(c*B).
Here a, b and c are the "parameters" of the model, B is named "the independent variable" and Y "the dependent variable".
Unfortunately you don't know what the values of a, b and c arre.
You then realize "observations" of Y for different values of B, and represent these latter by a vector named Bdata and the former by another vector of same length) named Hdata.
The procedure Fit is aimed to assess the values of a, b and c given the observations Bdata and Hdata.

Let Z(a, b, c) the vector whose component i is equal to a*Bdata[i]+b*sinh(c*Bdata[i]).
The most common method to seek for the triple (a, b, c) is to minimizes the euclidian norm of the vector Z(a, b, c) - Hdata wrt (a, b, c). This strategy is named "least-squares method".

I don't remember what Fit exactly did in Maple 13, but in more recent versions Fit is a procedure in the Statistics package which
"... Fit command fits a model function to data by minimizing the least-squares error...." (from help pages).

Your command has exactly the same syntax today:

the first parameter is the model,
the second one is the vector of the independent variable(s),
the third oneis the one of the dependent variable,
and the last one is the name of the independent variable (or the list of the independent variables).

Note that the command Fit can contain several options.

You write that you "never used maple13 and have to describe a code" does it mean you do not have Maple 13 right now?
If it is so, go here
https://www.maplesoft.com/documentation_center/history.aspx
select Maple 13 and next Maple User Manual to load the pdf.
Open it and search for fit or curve fitting.

No formal solution...

Answered: mmcdara 6149

April 19 2023

0 0

Use fsolve instead of solve, details in the attached file.

>

restart:

>	f := s -> cos(s)cosh(s); plot([1, f(s)], s=-4Pi..4*Pi, color=[red, blue], view=[default,-1/2..3/2])

>	# solve returns only the trivial solution s=0 and its not capable # get the (infinite) set of al solutions: infolevel[solve] := 4; solve(f(s)=1); infolevel[solve] := 0;

Main: Entering solver with 1 equation in 1 variable
Dispatch: dispatching to Trig handler
Dispatch: trigonometric substitution _S000002 = cos(s)
Recurse: recursively solving 1 equations in 2 variables
Dispatch: dispatching to OnlyIn handler
Recurse: recursively solving 1 equations in 2 variables
Dispatch: handling a single polynomial
Main: polynomial system split into 1 parts under preprocessing
Main: subsystem is essentially univariate
UnivariateHandler: subsystem has only one equation
UnivariateHandler: solving as if univariate in _S000002
UnivariateHandler: subsystem has only one equation
UnivariateHandler: solving as if univariate in _S000003
Polynomial: removing 1 trivial zeros at the origin
UnivariateHandler: subsystem has only one equation
UnivariateHandler: subsystem has only one equation
UnivariateHandler: solving as if univariate in _S000002
Polynomial: solving a linear polynomial
Recurse: recursively solving 1 equations in 1 variables
Dispatch: dispatching to OnlyIn handler
Recurse: recursively solving 1 equations in 1 variables
Transformer: solving the uncoupled linear subsystem in _S000004
Recurse: recursively solving 1 equations in 1 variables
Transformer: solving the uncoupled linear subsystem in s
SolutionsLost: setting solutions lost flag
Main: solving successful - now forming solutions
Main: Exiting solver returning 1 solution

Warning, solutions may have been lost

(1)

>	# f is symmetric: simplify(f(s)-f(-s))

(2)

>	# The roots of f are 0, Pi/2 + k*Pi where k is a relative integer roots_of_f := RealDomain:-solve(f(s)=0, allsolutions);

(3)

>	# In the half line [0, +oo), the first root is s=0 f(0)

(4)

>

# In the half line (Pi/2, +oo), f(s) > 0 when cos(s) > 0.
# This f(s) > in intervals (Pi/2+k*Pi, Pi/2+(k+1)*Pi) with k
# any odd integer.
# Let omega[k] such an interval.
# As cosh is a strictly increasing function over (0, +oo),
# there exist two points z1 and z2 omega[k] such that f(s) = 1.
#
# In the half line (-oo, 0), f(s) > 0 when cos(s) > 0.
# This f(s) > in intervals (Pi/2-k*Pi, Pi/2-(k+1)*Pi) with k
# any odd integer.
# Let omega[k] such an interval.
# As cosh is a strictly increasing function over (-oo, 0)),
# there exist two points z1 and z2 omega[k] such that f(s) = 1.

>

AllSolutionsInRange := proc(r)
  local kvalues, i, sols;
  op~({isolve({Pi/2+k*Pi >= max(0, op(1, r)), Pi/2+(k+1)*Pi <= max(0, op(2, r))})}):
  kvalues := rhs~(select((x -> rhs(x)::odd), %)):
  op~({isolve({Pi/2+k*Pi >= min(0, op(1, r)), Pi/2+(k+1)*Pi <= min(0, op(2, r))})}):
  kvalues := kvalues union rhs~(select((x -> rhs(x)::odd), %));
  sols := NULL;
  for i in kvalues do
    sols := sols, fsolve(f(s)=1, s=Pi/2+i*Pi..Pi/2+(i+1/2)*Pi);
    sols := sols, fsolve(f(s)=1, s=Pi/2+(i+1/2)*Pi..Pi/2+(i+1)*Pi);
  end do;
  if verify(0, r, interval) then sols := sols, 0 end if:
  return sort([sols])
end proc:

>	# Example r := -5Pi/2..9Pi/2; AllSolutionsInRange(r)

(5)

>

Download 1_mmcdara.mw

This is an illustration of what can be d...

Answered: mmcdara 6149

April 16 2023

2 2

This is an illustration of what can be done.
I used Maple 2015.2 but the structure which draw the label names is the same in bewer Maple's version.

Download Example.mw

An answer which complete my previous com...

Answered: mmcdara 6149

April 14 2023

1 0

The attached file presents:

A non linear model (NLM) and its linear counterpart (LM) under a simple transformation.
Two classical estimators of the covariance matrix from the NonlinearFit of NLM and how they both fail.
Two methods to asses the covariance matrix of the native parameters of the LinearFit of LM, from the one of the parameters of LM.
A discusion about the precautions to use as in both cases (NLM and LM) the joint distribution of the native parameters cannot be multi-normal.

Nonlinear_to_Linear_1.mw

Your initial question was "Need the covariance matrix...".
If you have read carefully the last text in the attached file you must have understood that in cas of a non linear model this covariance matrix doesn't provide sufficient information to draw inferences.
So my question is "Why do you need it?"

You probably need to add assumptions...

Answered: mmcdara 6149

April 13 2023

0 1

Example

Download total_PE_with_assumptions.mw

Don't know if it may help ...

Answered: mmcdara 6149

April 13 2023

3 6

Before several customization features have been introduced in the GraphTheory package, I used to do my own customization.
It can be lengthy but quite simple if you observe the structure that DrawGraph displays.

I don't know what Maple 2023 is capable to do, but in case it can do exactly what you want, there is an example of a method to change the shape, position and style of the arrows or the style of the vertices.
procedures MoveArrow and ChangeVertex have to be used after you have exploited all the possibilities offered by GraphTheory.

>

restart:

>	interface(version)

(1)

>	with(GraphTheory): with(plots):

>	G := Digraph({[1, 2], [2, 3], [3, 4], [4, 1]}):

>	dg := DrawGraph(G):

>	print~([op(dg)]):

(2)

CHANGING EDGES

>

# The last POLYGONS structures describe each edges with its arrow.
#
# If you want to change the style and position of any arrow you can do this.

MoveArrow := proc(edge, pc, angle, length, closed, col)
  local f, dx, dy, alpha, beta1, beta2, head, end1, end2, arrow, style:
  f      := (dx, dy) -> piecewise(dx >= 0, arctan(dy/dx), Pi+arctan(dy/dx)):
  dx, dy := (op(1, edge)[2] -~ op(1, edge)[1])[];
  alpha := f(dx, dy);
  beta1 := alpha-angle;
  beta2 := alpha+angle;
  head   := op(1, edge)[1]*~pc +~ op(1, edge)[2]*~(1-pc):
  end1   := head -~ length *~ [cos, sin](beta1):
  end2   := head -~ length *~ [cos, sin](beta2):
  style := remove(type, [op(edge)], list)[]:
  if closed then
    arrow := [end1, head, end2, end1]:
   #return POLYGONS([op(1, edge)[], ListTools:-Reverse(op(1, edge))[]], style),
   #       POLYGONS(arrow, COLOR(RGB, op(ColorTools:-NameToRGB24(col))))
    return POLYGONS(arrow, COLOR(RGB, op(ColorTools:-NameToRGB24(col)))),
           POLYGONS([op(1, edge)[], ListTools:-Reverse(op(1, edge))[]], style)
  else
    arrow := [end1, head], [head, end2]:
    return POLYGONS(op(1, edge), arrow, style)
  end if:
end proc:

NV     := numelems(Vertices(G)):
NE     := numelems(Edges(G)):
pol    := select(has, [op(dg)], POLYGONS):
edges := pol[-NE..-1]:
others := pol[1], remove(has, [op(dg)], POLYGONS)[]:

cols     :=["Red", "Cyan","Yellow", "Magenta"]:
pos      := [0.4, 0.3, 0.4, 0.3]:
NewEdges := seq(MoveArrow(edges[n], pos[n], Pi/12, 0.1, true, cols[n]), n=1..NE):

display(PLOT(others, NewEdges), axes=none)

CHANGING VERTICES

>	centers := map(n -> Statistics:-Mean(convert(op(n, others[1]), Matrix)), [$1..NV]): diameters := map(n -> abs~(op(n, others[1])[1] -~ centers[n]), [$1..NV]):

diameters := [Vector[row](2, {(1) = HFloat(0.0288), (2) = HFloat(0.028799999999999937)}), Vector[row](2, {(1) = HFloat(0.0288), (2) = HFloat(0.0288)}), Vector[row](2, {(1) = HFloat(0.028799999999999937), (2) = HFloat(0.028799999999999937)}), Vector[row](2, {(1) = HFloat(0.028799999999999937), (2) = HFloat(0.0288)})]

(3)

>

ChangeVertex := proc(shape, c, d, theta, col)
  local fig, n:
  if shape = "disk" then
    fig := display(plottools:-disk([0, 0], d[1], color=col));
  elif shape = "ellipse" then
    fig := display(plottools:-ellipse([0, 0], entries(d, nolist), color=col, filled=true))
  elif shape[1] = "P" then
    n   := parse(shape[2..-1]):
    fig := display(polygonplot([seq([cos(2*Pi*i/n), sin(2*Pi*i/n)], i = 1..n)], color=col)):
    fig := plottools:-scale(fig, entries(d, nolist))
  end if;
  plottools:-translate(plottools:-rotate(fig, theta), entries(c, nolist));
end proc:

>

# shapes := {"disk", "ellipse", ...) (see plottools)
#           or "Pn" where n is an integer >= 3

NewDiameters := map(d -> d*~[2, 1], diameters):
rotations    := (Pi/NV)*~[$1..NV]:
shapes       := ["disk", "ellipse", "P4", "P7"]:
display(
  seq(ChangeVertex(shapes[n], centers[n], NewDiameters[n], rotations[n], "Chartreuse"), n=1..NV)
  , PLOT(others[2..-1], NewEdges)
  , axes=none
  , scaling=constrained
)

>

Download Customization_1.mw

I don't have time right now to answe...

Answered: mmcdara 6149

April 13 2023

1 5

I don't have time right now to answer your second question.
More of this it requires some precisions: for instance what are S[1] and S[2] for they do notappear in the expresison of Omega?

First question:
You have made an error in reasoning: see the explanation in the attached file after the yellow highlighted text.
Basically the variance and covariance are sufficient statistics from which other statistics can be build (the linear correlation coefficient for instance), not the opposite.
Thus you cannot do some transformation on random variables, believe the correlation coefficient is an invariant, and deduce the covariance from its value and the two variances.

RVs_sum_mmcdara.mw

About sufficient statistics see here https://en.wikipedia.org/wiki/Sufficient_statistic

Rescaling is not the very problem...

Answered: mmcdara 6149

April 11 2023

2 0

I agree with @dharr that it isalways better to rescale the data to avoid numerical problem.
But rescaling hides a major bug in Fit/LinearFit: both these functions say C1=0 while it's not.

I already observed several weaknesses or failures of LinearFit, this is just a new one.
In case you doubt the results LinearFit returns, use the formula to get the vector A of the regressor coefficients:

A := (Z^+ . Z)^(-1) . Z^+ . Y;

In the attached file you will see that the wrong result returned by Fit/LinearFit comes from a wrong estimation of the rnk of the information matrix ((Z^+ . Z) in the text): as this rank is rqual to 3 Fit/LinearFit finds it equal to 2 only for precision problems.

These issues can be solve by rescaling the data,; but rescaling them BECAUSE C1=0 is not a good argument.

>

restart:

>	with(Statistics):

>	interface(version)

(1)

>

Data := Matrix([[4.74593554708566, 11385427.62, 2735660038000], [4.58252830679671, 25469809.77, 12833885700000], [4.29311160501838, 1079325200, 11411813200000000], [4.24176959154225, 1428647556, 18918585950000000], [5.17263072694618, 1428647556, 18918585950000000], [4.39351114955735, 1877950416, 30746202150000000], [4.39599006758777, 1428647556, 18918585950000000], [5.79317412396815, 2448320309, 49065217290000000], [4.48293612651735, 2448320309, 49065217290000000], [4.19990181982522, 2448320309, 49065217290000000], [5.73518217699046, 1856333905, 30648714900000000], [4.67943831980476, 3071210420, 75995866910000000], [4.215240105336, 2390089264, 48670072110000000], [4.41566877563247, 3049877383, 75854074610000000], [4.77780395369828, 2910469403, 74061327950000000], [4.96617430604669, 1416936352, 18891734280000000], [4.36131111330988, 1416936352, 18891734280000000], [5.17783192063198, 1079325200, 11411813200000000], [4.998266287191, 1067513353, 11402362980000000], [4.23366152474871, 2389517120, 48661380410000000], [4.58252830679671, 758079709.3, 5636151969000000], [6.82390874094432, 1304393838, 14240754750000000], [4.24176959154225, 912963601.2, 8621914602000000], [4.52432881167557, 573965555.4, 3535351888000000], [4.84133601918601, 573965555.4, 3535351888000000], [6.88605664769316, 732571773.2, 5558875538000000], [4.35575841415627, 1203944381, 13430693320000000], [4.42527441640593, 955277678, 8795128298000000], [6.82390874094432, 997591754.9, 8968341995000000], [4.35144484433733, 143039477.1, 305355143300000]]):

>	X := Data[.., [2, 3]]: Y := Data[.., 1]: Fit(C1+C2v+C3w, X, Y, [v, w]); LinearFit(C1+C2v+C3w, X, Y, [v, w]);

Warning, model is not of full rank

(2)

>	Z := `<\|>`(Vector(numelems(Y), 1), X): r := proc(d) Digits:=d: printf("Digits %2d rank = %d\n", d, LinearAlgebra:-Rank(Z^+ . Z)); Digits:=10 end proc: r(10), r(20), r(40): Digits(); A := (Z^+ . Z)^(-1) . Z^+ . Y; F := `<,>`(1, x1, x2) . A assuming x1::real, x2::real;

Digits 10 rank = 1
Digits 20 rank = 2
Digits 40 rank = 3

A := Vector(3, {(1) = 4.659353816079307, (2) = 0.5985084089529947e-9, (3) = -0.27350964718426345e-16}, datatype = float[8])

(3)

>	plots:-display( ScatterPlot3D(`<\|>`(X, Y), symbol=solidsphere, color=blue, symbolsize=20), plot3d(F, x1=(min..max)(X[..,1]), x2=(min..max)(X[..,2]), style=wireframe, color=gray) );

WITH SCALED DATA

>	model := unapply( LinearFit(C1+C2v+C3w, Scale(X), Scale(Y), [v, w]), [v, w])

(4)

>	mX, sX := (Mean, StandardDeviation)(X); mY, sY := (Mean, StandardDeviation)(Y);

(5)

>	MODEL := model((x1-mX[1])/sX[1], (x2-mX[2])/sX[2]) * sY + mY

(6)

>	plots:-display( ScatterPlot3D(`<\|>`(X, Y), symbol=solidsphere, color=blue, symbolsize=20), plot3d(F, x1=(min..max)(X[..,1]), x2=(min..max)(X[..,2]), style=wireframe, color=gray), plot3d(MODEL, x1=(min..max)(X[..,1]), x2=(min..max)(X[..,2]), style=wireframe, color=red) );

>

Download Bug_In_Fit.mw

The explicit computation of A and LinearFit applied to Scaled data give the same result.
At the minimum the fact that LinearFit applied on raw data gives a wrong result is a weakness, but I lean towards a bug:

either LinearFit must automatically rescale the data,
eiither LinearFit must say "The information matrix is numerically not of full rank, so I don't return any result".

Note: the only possibility for this matrix not to be of full rank (3) is that there exist three numbers a, b, c such that.

a*X[..,1]+b*X[..,2]+c = 0

It's easy to verify this is not true:

ScatterPlot(X[..,1], X[..,2]);
# the residual sum of squares = 0 iif a*x1+b*x2+c = 0
LinearFit([1, x1], X[..,1], X[..,2], [x1], output=residualsumofsquares)
                     1.03756169900198*10^33
# other way
LinearAlgebra:-LinearSolve(Z, Vector(numelems(Y), 0))
               < a, b, c > =  < 0, 0, 0 >

It's simpler and less sensitive to roundoff errors to use this later command to check if the matrix (Z^+ . Z) is of full rank, or nor.
Note that you can use identically

LinearAlgebra:-LinearSolve(X, Vector(numelems(Y), 0))

At least compare the results these 4 commands return

LinearAlgebra:-Rank(Z^+ . Z);         # uncorrect result (should be 3)
LinearAlgebra:-Rank(Z);               # uncorrect result (should be 3)
                               1
                               2
LinearAlgebra:-Rank(X^+ . X);         # uncorrect result (should be 1)
LinearAlgebra:-Rank(X);               # correct result
                               1
                               2

This does not add much to what has been ...

Answered: mmcdara 6149

April 11 2023

3 0

This does not add much to what has been said before.
I just imagined that the challenge is to find "by hand" a formal development of P_n(x) ????

So the idea I followed is this one:

Start from n = 1 and interpret your equality this way: P₁(x) is the quotient of (x-x^₂)⁴ by (1+x²) and -4 is its remainder.
The division of (x-x^₂)⁴ by (1+x²) can easily be done by hand to find the expression of P₁(x) and the remainder -4.
The remainder of the division of (x-x^₂)⁴ by (1+x²) is obviously (-4)ⁿ .
The binomial expansion then gives the expression of P_n(x) as a function of P₁(x) and r₁=-4.

>

restart

>	f := n -> (x-x^2)^(4*n); g := (1+x^2);

(1)

>	# You have this relation f(n) = g*p[n] + (-4)^n

(2)

>	# An interpretation is "the ramainder of the division of f(n) by # (1+x^2) is (-4)^n and its quptient is (1-x)*p(n)(x): # # For n = 1: do the euclidian division of f(1) by (1+x^2) QR := proc(n) quo(f(n), (1+x^2), x, 'r'): return %, r end proc: p[1], r[1] := QR(1)

(3)

>	# then 'f(1)' = 'p[1]'*(1+x^2) + r[1]

(4)

>	# For n > 1 f(n) writestat 'f(n)' = ('p[1]'*(1+x^2) + 'r[1]')^n

(5)

>	# Thus the remainder of division of f(n) by (1+x^2) # is r[n] := r[1]^n

(6)

>	# And its quotient q[n] verifies ('p[1]'(1+x^2) + 'r[1]')^n = (1+x^2)p[n] + 'r[1]'^n

(7)

>	# Using the binomial expansion of the lhs from k=1 to k=n # (to remove the constant term) one gets the expression of # q[n] = P(n): P := n -> sum(factorial(n)/factorial(k)/factorial(n-k)(1+x^2)^(k-1)p[1]^k*r[1]^(n-k), k=1..n)

proc (n) options operator, arrow; sum(factorial(n)*(x^2+1)^(k-1)*p[1]^k*r[1]^(n-k)/(factorial(k)*factorial(n-k)), k = 1 .. n) end proc

(8)

>

P(1);

(9)

>	expand(P(2))

(10)

>

Download By_hand.mw

A way...

Answered: mmcdara 6149

April 09 2023

2 1

This is not exactly what you want, but it's not that far.
You can still tune this code to improce the rendering, in particulat the positions and fonts of the three parameter counters.

>

restart;

>	with(DEtools):with(plots):with(plottools):

>	sigma1:=e1alpha: sigma2:=e2delta:

>	g:=x/(1+betax^2); f:=(thetax-1)(1-x)(1+betax^2)-y; h:=alphag-z-sigma1; j:=delta*y-sigma2;

alpha*x/(beta*x^2+1)-e1*alpha-z

(1)

>

anim := proc()
  local p0, p1, p3, q0, q1, q2, q3:

  p0:=theta->display(
               plot([1/theta,y,y=0..1],linestyle=dash,color=green),
               textplot([0, 1.05, `#mo("θ")`=evalf[3](theta)], align=right)
             ):
  p1:=e1->display(
               plot([x,e1,x=0....1.5],color=blue),
               textplot([1.5, 1.05, `#mo("e1")`=evalf[3](e1)], align=left)
             ):
  p3:=beta->display(
               plot([x,x/(1+beta*x^2),x=0..1.5],color=magenta),
               textplot([0.75, 1.05, `#mo("β")`=evalf[3](beta)])
             ):

  q0:=animate(p0,[theta],theta=2...10):
  q1:=animate(p1,[e1],e1=0.1..1):
  q3:=animate(p3,[beta],beta=0..1.5):
  q2:=plot([1,y,y=0..1],linestyle=dash,color=green):

  [q0, q1, q2, q3]
end proc:

display(anim(),title=" ", view=[0..1.5,0..1.1])

>

Download animate_mmcdara.mw

For fun...maybe?...

Answered: mmcdara 6149

April 09 2023

0 3

Not cuneiform symbols but Tibetan marks.

Download For_fun.mw

Here is a list of some special ASCII characters: http://www.addressmunger.com/special_ascii_characters/
Maple can visualize (almost) all of these characters up to

𝐀

Capital A bold; Capital 𝐀 bold

𝐀

excluded

plots:-intersectplot...

Answered: mmcdara 6149

April 09 2023

1 0

But this requires to build the cartesian equation of the two tori:

with(plottools):
with(plots):

# torus([1,1,1], 1, 2)
# torus([1,6,1], 1, 2)

# Torus with symmetry plane // (O,x, yj) 
#   center       = [p, q, r], 
#   major radius = a, 
#   minor radius = b.
#
# 
#
# Cartesian equation

ce := (p, q, r, a, b) ->  ((x-p)^2 + (y-q)^2 + (z-r)^2 + a^2 - b^2)^2 = 4*a^2*((x-p)^2 + (y-q)^2):
  
# Cartesian equations of torus 1 and torus 2
T1 := ce(1, 1, 1, 2, 1);
T2 := ce(1, 6, 1, 2, 1);

# Their intersection is contained in the plane y=7/2

plots:-intersectplot(T1, T2, x=-2..4, y=3.4..3.6, z=0..2, color=red, thickness=5, grid=[20, 2, 20])

Download Intersection.mw

restart:with(Statistics):X := Vector([1,...

Answered: mmcdara 6149

April 09 2023

2 0

restart:
with(Statistics):
X := Vector([1, 2, 3, 4, 5, 6], datatype=float):
Y := Vector([2, 5.6, 8.2, 20.5, 40.0, 95.0], datatype=float):
f := unapply(PowerFit(X, Y, v), v)
                                 HFloat(2.041674366793877)
v -> HFloat(1.4795777879958958) v                         
f(Pi)
                   HFloat(15.316373717844828)

Why do you use a deprecated function?Not...

Answered: mmcdara 6149

April 08 2023

0 3

Why do you use a deprecated function?
Note that using solve instead doesn't change things, but at least you avoid unnecessary problems.

The fact is that the right and left limit of your expression at x=1/2 are not the same. More precisely they real parts are thesame, but they have opposite imaginary parts.

>

restart:

>	interface(version)

(1)

>	expr := (x-1/2)(1-sqrt(1+3x(x-2/3)/(5(x-1/2)^2)))/(x-2/3): sol := solve(expr);

(2)

>	limit(expr, x=sol[1], left); limit(expr, x=sol[1], right); convert(series(expr, x=1/2, 2), polynom);

((3/5)*I)*5^(1/2)*csgn(I/(x-1/2))+(-6+((12/5)*I)*5^(1/2)*csgn(I/(x-1/2)))*(x-1/2)

(3)

>	limit(expr, x=sol[2], left); limit(expr, x=sol[2], right);

(4)

>	plots:-dualaxisplot( plot(Re(expr),x=-1/4..3/4, color=blue, numpoints=400, thickness=3, discont=true, legend="Re"), plot(Im(expr),x=-1/4..3/4, color=red, numpoints=400, thickness=3, discont=true, legend="Im"), gridlines=true )

>

Download solve_Re_Im.mw

E-Mail Address:
Password:
Remember Me:	Automatically sign in on future visits

E-Mail Address:
Password:
Remember Me:	Automatically sign in on future visits

Ask a Question

Create a Post

6149 Reputation

17 Badges

MaplePrimes Activity

These are answers submitted by mmcdara

A question ofsize....

You have two different data, let's s...

No formal solution...

This is an illustration of what can be d...

An answer which complete my previous com...

You probably need to add assumptions...

Don't know if it may help ...

I don't have time right now to answe...

Rescaling is not the very problem...

This does not add much to what has been ...

A way...

For fun...maybe?...

plots:-intersectplot...

restart:with(Statistics):X := Vector([1,...

Why do you use a deprecated function?Not...

Save this setting as your default sorting preference?

Ask a Question

Create a Post

Generating PDF…

Save this setting as your default sorting preference?
Note: You can change your preference any time in your account settings
Don't show this again