a:=ifactor(op(1, convert(0.999987406876435, fraction)));

b:=ifactor(op(2, convert(0.999987406876435, fraction)));

c:=ifactor(op(1, convert(0.999919848203811, fraction)));

d:=ifactor(op(2, convert(0.999919848203811, fraction))); 

x1:=2:

x2:=3:

x3:=7:

x4:=17:

x5:=173:

x6:=709:

x7:=5347:

x8:=18713:

L:=[seq([x||i, x[i]], i=1..8)]:

for k in [a, b, c, d] do

subs({seq(L[i,1]=cat('x', op(L[i,2])), i=1..8)}, k);

od;

convert / phaseamp  appeared only in recent versions of Maple. If you have an older version, you can use  applyrule  command:

Rule:=a::realcons*sin(x)+b::realcons*cos(x)=sqrt(a^2+b^2)*sin(x+arctan(b,a)):

applyrule(Rule, sin(x)+cos(x));

To plot a graph of  phi  against time is  easier than making of animation. This can be done directly by the command  plots:-odeplot  using the procedure  integ :

ode := {diff(phi(t), t, t) = -phi(t), phi(0) = (1/12)*Pi, (D(phi))(0) = 0};

integ := dsolve(ode, numeric);

plots:-odeplot(integ, [t, phi(t)], 0..15);

We assume that  A - is a known matrix, for example,  A = <<1 | -2 | 3> , <4 | 1 | 2> , <5 | 1 | 2>> ,  and  B - is a subspace of  R^3, for example, all vectors collinear to  <1, 2, 3> . First, we find the general form of the matrices  Y  for which  NullSpace(Y) = B :

restart;

Y:=<<a[1] | b[1] | c[1]> , <a[2] | b[2] | c[2]> , <a[3] | b[3] | c[3]>>:

solve({op(Equate(Y.<1, 2, 3>, <0, 0, 0>))}, {a[1], a[2], a[3]});  assign(%);

Y;

Columns  <b[1], b[2], b[3]>  and  <c[1], c[2], c[3]>  should not be proportional.

Next, we find the matrix  X :

X:=Matrix(3, symbol=x);

A:=<<1 | -2 | 3> , <4 | 1 | 2> , <5 | 1 | 2>>;

Sys:=op(Equate(A.X+X.A, Y));

solve({Sys}, {seq(seq(x[i,j], j=1..3), i=1..3)});

We see that there are infinitely many solutions, depending on the 6 parameters.

You can use  RealRange command . For example:

is(2 in RealRange(0, Open(2)));

is(1 in RealRange(0, Open(2)));

is(0 in RealRange(0, Open(2)));

                              false

                              true

                              true

When solving equations, you can just use the corresponding inequalities. Because in your equation 3 unknowns, two of them must be specified. For example:

restart;

y:=2:  z:=0.3:

solve({y=1.048 + 1.02*x + 6.118*(z-4.041*x^2) + 16.22*(z^2) +6.241* (x*z),

x>=0.001, x<=0.7}, x);

                                                         {x =0 .3718638796}

a := x->piecewise(x^2 <= 3, x^2, x^1.5):

a(1);  a(2);

                  1

        2.828427125

Addition: one more way

a := x->`if`(x^2 <= 3, x^2, x^1.5):

convert(exp(x), Sum, dummy=n);

This is quite easy. For example:

V:=[V[k] $ k=1..n]: N:=[$ 1..n]:

plot(N, V, thickness=2);

Even for a specified  b the  solution is expressed only in terms of integrals:

for b from -5 to 5 do

P[b](t):=t*ln(t)^b:

dsolve(diff(Q(t), t)^2 = diff(P[b](t), t), Q(t)):

od;

1st column of  matrix  A  by  3rd column of matrix  B :

A := Matrix([[1, 2, 3], [4, 5, 6], [7, 8, 9]]);
B := Matrix([[9, 8, 7], [6, 5, 4], [3, 2, 1]]);

A[.., 1] .~ B[.., 3];  whattype(%);

Edited: the extra bracket was deleted. 

restart;

with(LinearAlgebra):

A := `<|>`(`<,>`(1, 2,3), `<,>`(2, 3, 0), `<,>`(2, 0, 0));

v, EigenVector1:= Eigenvectors(A):

NewMatrix3 := Matrix([[x1,x2,x3], [x4,x5,x6], [x7,x8,x9]]):

simplify(solve({seq(seq((NewMatrix3 . EigenVector1(1..-1,i))[j]=(v[i]*EigenVector1(1..-1,i))[j], j=1..3), i=1..3)}, {seq(x||i, i=1..9)}));

assign(%);

NewMatrix3;

 Edited:  restart  command was added.

L := [3, -2, 5, 1, 10, 7]:

[seq(add(L[i], i = 1 .. k), k = 1 .. nops(L))];

                   3, 1, 6, 7, 17, 24

Addition:  if  L  is a long list, then more effective variant is

L := [3, -2, 5, 1, 10, 7]:  S[1] := L[1]:  n := nops(L):

for i from 2 to n do  S[i] := S[i-1]+L[i]  end do:

S := [seq(S[i], i = 1 .. n)];

I just built your plot by points (a total of 41 points were used):

f:=phi->Int(309.9615685*exp(3.936759889*phi-3.936759889*`&varphi;`)*Heaviside(-phi+`&varphi;`)*exp(-1549.807842*`&varphi;`^2)/(Pi*sqrt(-15.74703956*phi+15.74703956*`&varphi;`)), `&varphi;` = -infinity .. infinity);

X:=[seq(x, x=-1..1, 0.05)]:

Y:=[seq(f(X[i]), i=1..nops(X))]:

plot(X, Y, thickness=2);

assume(L > 0);
assume(K > 0);
assume(mu > 0);
assume(a<0);
assume(t > 0);
convert(cosh((1/2)*t*sqrt(a)/L), sincos);
subs(a=mu^2-4*L*K, %);

Code is edited.