# slelievre's profile - activity

2021-01-26 02:17:42 -0600 edited answer Solving the Schrodinger equation using sage (Well)

Hello, @Grim! I am not quite sure about Schrödinger's Equation, but I know how to obtain $k=\frac{n\pi}{a}$. I believe this will result a helpful version of your code:

k, a = var('k, a', domain='positive')
psi = function('psi')(x)
eqs = diff(psi, x, 2) + k^(2)*psi == 0
EQS = desolve(eqs, psi, ivar=x)
A, B = var('A, B')
EQS = EQS.subs(_K1=B, _K2=A)
show(LatexExpr(r'\psi(x) = '), EQS)
eq1 = EQS.subs(x=0) == 0  # this is psi(0) = 0
eq2 = EQS.subs(x=a) == 0  # this is psi(a) = 0
show(eq1)
show(eq2)


The first line defines the symbolic variables k and a to be positive, so you can avoid using the assume command. The second line defines psi as the unknown function; the third line defines the ordinary differential equation; and the fourth one solves the differential equation.

Now, things get a little complicated. If you write print(EQS), you will notice that there are a couple of new variables called _K1 and _K2. They are integration variables (byproducts of the resolution of the ODE). However, they are not symbolic variables that you have defined, so you won't be able to use them outside of EQS (for example, to use the solve command or any other). So, I replaced them with A and B in the sixth line. The seventh line (the first show command) is there only to check that I have obtained the same definition of $\psi$ as you in your code (it will show you $\psi(x)=A\cos(kx)+B\sin(kx)$). Then, eq1 and eq2 represent your boundary conditions $\psi(0)=0$ and $\psi(a)=0$, respectively.

At this point, there are many ways in which you can obtain the values of $A$, $B$ and $k$. The first method is manually solving equations eq1 and eq2, which are shown by the last two lines of my code. You will see that

$$A=0$$ and $$A\cos(ak)+B\sin(ak)=0$$

If you replace the first equation in the second, you obtain $B\sin(ax)=0$. Now, since $B>0$ (according to your code), you must have that $\sin(ak)=0$, which implies that $k=\frac{n\pi}{a}$. Finally, in this case, $B$ can be any positive number (regardless of its value, your ODE and boundary conditions are satisfied.)

The second way to solve this is using Sage itself. You have two equations, namely eq1 and eq2, so you will need two variables. You can write

sol = solve([eq1, eq2], A, k)
show(sol)


just after the previous code. You will get something like

$$\left[\left[A = 0, k = \frac{\pi + 2\pi z_{1879}}{a}\right], \left[A = 0, k = \frac{2 \pi z_{1900}}{a}\right]\right]$$

The variables z1879 and z1900 are similar to the _K1 and _K2 that I mentioned above. More precisely, they are integer parameters; just mentally replace them with $m$ and $n$. Then you have two possible solutions. After further analysis, you will notice that these two solutions can be written more simply as $A=0$ and $k=\frac{n\pi}{a}$. Once again, $B$ has no restrictions, besides being positive.

You can also use the Sympy package to solve these equations:

sol = solve([eq1, eq2], A, k, algorithm='sympy')
show(sol)


You will get

$$\left[\left\{A : 0, k : 0\right\}, \left\{A : 0, k : \frac{\pi}{a}\right\}\right]$$

These are just two particular solutions that you can easily generalize using reduction formulas and properties of the trigonometric functions.

What happens if you use solve([eq1, eq2], A, B) or solve([eq1, eq2], B, k)? You will get $[[A=0, B=0]]$ and $[]$, respectively. The first is a correct solution to both equations, but not one that you are looking for. The second is Sage telling you that could not find a solution. If you use SymPy, you will get $\left[\left\{B: -\frac{A}{\tan(ak)}\right\}\right]$, which is a way of admission of defeat, but it is more useful, because it hints you that it is Sage's (actually SymPy's) fault, while the $[]$ before could lead you to think that there is no solution.

I would like to point out a few interesting teachings derived from your question:

1. Do not depend on any software to solve a problem; it's just a tool.
2. SageMath is excellent in solving problems, but it is not perfect. Try different approaches (e.g., as we did when we used SymPy), especially if the problem is complex.
3. You can use a software to help you, but the more of the resolution you can do manually (using ingenuity), the more sure you will be of your result. In the end, if you can do it manually, why bother using Sage? In the end, there are things that are easy for a human mind, but very difficult for a software.
2021-01-26 02:10:30 -0600 commented answer Solving the Schrodinger equation using sage (Well)

@dsejas -- this might be due to a bug in the way LaTeX is processed on Ask Sage.

It seems some of the backslashes need escaping, ie we need \\ instead of \.

In my experience, not all backslashes need this, but if you do them all you're safe.

It seems to depend on the following character, and seems

• required for \{, \}, \p
• not required for \c, \s, \t

and one could of course prepare a more complete list.

The preview when composing a question or answer behaves differently, making this all the more puzzling.

2021-01-25 13:53:24 -0600 edited question Solving the Schrodinger equation using sage (Well)

I'm trying to solve the Schrödinger equation through Sage, for the potential well.

I wrote:

k, a = var('k, a')
psi = function('psi')(x)
eqs = diff(psi, x, 2) + k^(2)*psi == 0
assume(k > 0, a > 0)
EQS = function('EQS')(x)
EQS = desolve(eqs, psi, ivar=x)
A, B, c = var('A, B, c')
result0 = A*cos(k*0) + B*sin(k*0) == 0
result1 = B*sin(k*a) == 0
assume(B > 0)


How can i find the constants A and B through boundary conditions?

My conditions are: Psi(0)=0 and Psi(a)=0.

If i put it in desolve function it return me zero.

So i tried to work with the equations separated but i cannot figure out how to get k=n*pi/a (the expected result).

2021-01-25 13:50:19 -0600 commented question Trying to solve the Schrodinger equation using sage (well)
2021-01-25 13:31:55 -0600 answered a question Arithmetic operation is not working in finite field GF(4091^2)?

The difference comes from different methods being available depending on the underlying finite field implementation.

The givaro implementation, available only for prime powers less than 2^16, provides a fetch_int method.

Below we illustrate this and propose a fetch_int function to provide the functionality in other cases.

Illustration

sage: p = 13
sage: F.<a> = GF(p^2)
sage: b = F.fetch_int(150)
sage: c = F.fetch_int(97)
sage: b, c, b + c
(11*a + 7, 7*a + 6, 5*a)

sage: pp = 4091
sage: FF.<aa> = GF(pp^2)
sage: FF.fetch_int(150)

sage: type(F)
<class 'sage.rings.finite_rings.finite_field_givaro.FiniteField_givaro_with_category'>
sage: type(FF)
<class 'sage.rings.finite_rings.finite_field_pari_ffelt.FiniteField_pari_ffelt_with_category'>

sage: FF.<aa> = GF(pp^2, impl='givaro')
Traceback (most recent call last)
...
ValueError: q must be < 2^16

sage: F.fetch_int??
...
Given an integer n return a finite field element in self
which equals n under the condition that :meth:gen() is set to
:meth:characteristic().
...


Define fetch_int as a function:

def fetch_int(F, n):
return FF(n.digits(base=F.characteristic()))


Then use it as follows:

sage: fetch_int(FF, 12275)
3*a + 2
sage: 3*pp + 2
sage: n
12275


Making all finite fields provide fetch_int regardless of the implementation is the object of

2021-01-24 16:44:38 -0600 edited answer Generating only the ranks of elliptic curves that can be found provably correctly

Here is a piece of code trying to compute the rank of all elliptic curves $E(a,b)$ of the shape $y^2=x^3+ax+b$ for $a,b\in[2010, 2021]$. Sometimes, there it will be "harder" to compute the rank, so the code will not deliver a computed rank. We fill in a dictionary with keys $(a,b)$ and values the corresponding rank for the key when it could be computed, and None otherwise.

R, dic = [2010 .. 2021], {}
for a, b in cartesian_product([R, R]):
try:
E = EllipticCurve(QQ, [a, b])
E.two_descent(second_limit=13, verbose=False)
r = E.rank(only_use_mwrank=False)
dic[(a, b)] = r
print(f'({a}, {b}) -> {r}')
except Exception:
dic[(a, b)] = None


To see which cases could not be computed...

[key for key, val in dic.items() if val is None]


... and get an empty list, so in this case there was no None value. But well, in other cases we may have the situation. (Also maybe consider first commenting out the two_descent line, as it may be time consuming; but it may allow finding more curve ranks.)

2021-01-22 11:37:49 -0600 commented answer How to do low degree computation in a Free Algebra ?

2021-01-22 11:30:42 -0600 edited answer How to do low degree computation in a Free Algebra ?

Another way (less complete than @slelievre's answer but more efficient for just a quick computation) is to compute directly in the free algebra, using the following instead of the product:

def fast_product(a, b):
res = F.zero()
a = F(a)
b = F(b)
data_a = [(w.to_word(), cf) for w, cf in a]
data_b = [(w.to_word(), cf) for w, cf in b]
for wa, cfa  in data_a:
for wb, cfb in data_b:
if len(wb) + len(wa) <= k:
res += cfa * cfb * F.monomial(wa) * F.monomial(wb)
return res

2021-01-21 18:36:25 -0600 edited answer Generating only the ranks of elliptic curves that can be found provably correctly

I assume that you mean elliptic curves over ℚ.

I would do two or three things. First, wrap your call to E.rank() in a try/except block to catch cases where an error is raised. Second, create an E.mwrank_curve() object, which you can ask its rank and ask whether the result is proved correct (this is delivered by mwrank):

sage: E = EllipticCurve([0,0,1,-7,6])
sage: Em = E.mwrank_curve()
sage: Em.rank()
3
sage: Em.certain()
True

sage: E = EllipticCurve([0, -1, 1, -929, -10595])
sage: Em = E.mwrank_curve()
sage: Em.rank()
0
sage: Em.certain()
False


In the second example, the rank really is 0 but E has Sha of order 4 and mwrank is not able to tell the difference.

I hope this helps.

2021-01-21 16:46:14 -0600 edited answer How do I understand the result of symbolic integrals

Maxima 5.21.1 gives

-1/x-gamma\_incomplete(0,-x)-gamma\_incomplete(-1,-x)


for

integrate(diff((exp(x) - 1)/x, x), x)


which is correct from what I can tell (agrees numerically with the original expression and has the same derivative).

The result isn't as simple as it could be because the integration algorithm is phrased in more general terms, such that the integrand you specified is a special case of some general form. Often that's the most effective way to calculate integrals, since you can cover a lot of special cases with one general form.

2021-01-21 11:33:57 -0600 answered a question Mysterious "must have a value" error

The command you want requires the GAP package HAPprime, no longer being developed: its repository at GitHub is archived, and some of the code has been migrated to HAP.

That said, you can still install it by downloading the package's compressed tarball and unpacking it in your ~/.gap/pkg/ folder (after creating it if necessary).

To do that on CoCalc, open a CoCalc terminal and run the following commands (skip the initial $ on each line which represents the prompt): $ GAP_PKG='~/.gap/pkg'
$mkdir -p$GAP_PKG
$cd$GAP_PKG
$URL=https://files.gap-system.org/gap48/tar.gz/packages$ TGZ=happrime-0.6.tar.gz
$wget$URL/$TGZ$ tar xf $TGZ  Then you can use the desired GAP command, either directly in GAP or in Sage (via gap or, better, libgap) after loading the package HAPprime. In Sage: sage: libgap.LoadPackage("happrime") #I HAP warning: Set BROWSER_PATH manually if needed. true sage: G = libgap.SmallGroup(64, 134) sage: R = libgap.ResolutionPrimePowerGroupRadical(G, 10) sage: R Resolution of length 10 in characteristic 2 for <pc group of size 64 with 6 generators> . No contracting homotopy available. A partial contracting homotopy is available.  2021-01-21 05:01:59 -0600 commented answer Issues with Reed-Solomon encoder Note about zero and one: • Fx.zero() directly provides the zero element of Fx while Fx(0) uses the integer zero and then turns it into an element of Fx. • Likewise, Fx.one() directly provides the one element of Fx while Fx(1) turns the integer one into an element of Fx. So one can use E.encode(Fx.one()) and E.encode(Fx.zero()). One could name the polynomials 0, 1 and a for easy access, e.g. p0 = Fx.zero() p1 = Fx.one() pa = Fx(a)  and then use e.g. E.encode(pa) instead of E.encode(Fx(a)). 2021-01-20 12:18:55 -0600 commented question Solve returns bad results Welcome to Ask Sage! Thank you for your question! 2021-01-20 12:17:13 -0600 commented question Transforming a list of list Please provide self-contained code that can be copied and pasted. Currently, we get name errors, as neither A, B, C, Rien are defined. Calling A the table or matrix when one of the duelists is also A seems error-prone. 2021-01-20 12:03:42 -0600 commented question Tensor product of elements of non-free algebras Thanks for reporting. This appears to reveal a bug. Fixing it is now tracked at 2021-01-20 12:02:23 -0600 answered a question Tensor product of elements of non-free algebras Thanks for reporting. This appears to be a bug. Fixing it is now tracked at Note that the code a = SteenrodAlgebra(2).an_element() M = CombinatorialFreeModule(GF(2), ['s', 't', 'u']) s = M.basis()['s'] T = tensor([a, s])  works fine in Sage 8.8 but fails in Sage 8.9. This change might have to do with the changes in Sage Trac ticket 25603. 2021-01-20 07:47:12 -0600 commented question Extracting coefficients Indentation is fixed now. 2021-01-20 07:46:44 -0600 edited question Extracting coefficients I have defined and generated a system like below: def p(*t): return SR.var(('p' + '_{}' * len(t)).format(*t)) f = [[[sum(p(i, a, k, b) for b in (1 .. 2)) for a in (1 .. 2)] for k in (1 .. 2)] for i in (1 .. 2)] b = [p(a, b) for a in (1 .. 2) for b in (1 .. 2) for k in (1 .. 2)] v = vector([f[i][j][k] for i in (0 .. 1) for j in (0 .. 1) for k in (0 .. 1)]) e = [bi - vi ==0 for bi, vi in zip(b, v)]  I would like to know if it is possible to extract the matrix coefficients of the equations in 'e'. The difficulty comes from the definition of the variables 'p()'. For example: in 'e', I have 8 expressions as below: e: [p_1_1 - p_1_1_1_1==0, p_1_1 - p_1_2_1_1==0, p_1_2 - p_1_1_2_1==0, p_1_2 - p_1_2_2_1==0, p_2_1 - p_2_1_1_1==0, p_2_1 - p_2_2_1_1==0, p_2_2 - p_2_1_2_1==0, p_2_2 - p_2_2_2_1==0]  So, my question is if it is possible to extract the coefficients of all the variables p_ _ and p_ _ _ _ considering the lexicography order and taking the p's with 2 entrees as the firsts variables, i.e., with the order: p_1_1, p_1_2, p_2_1, p_2_2, p_1_1_1_1, p_1_1_1_2, ... and so on. It means, that I have a system with 4 variables with 2 entrees p_ _ and more 16 with 4 entrees p_ _ _ _. I would like to write these coefficients using matrices in this case with 20 variables and 8 equations given in e. Precisely, I want to know the degrees of freedom of the system. 2021-01-20 07:30:28 -0600 answered a question Transforming a list of list The function individual_duel_scores below returns a dictionary of dictionaries with the duel scores. One can use it to build a matrix for all the duel scores. The function: def individual_duel_scores(duel_winners, duel_scores): r""" Return a dictionary of dictionaries giving the scores of each players against each other player. """ from collections import defaultdict players = sorted(set(flatten(duel_winners))) scores = defaultdict(lambda: defaultdict(int)) for ((a, b), _), (x, y) in zip(duel_winners, duel_scores): scores[a][b] = x scores[b][a] = y return scores  Define the players: sage: A, B, C = 'ABC'  The duel winners and the scores (called l1 and l2 in the question): sage: duel_winners = [([𝙰, 𝙱], 𝙰), ([𝙰, 𝙲], 𝙲), ([𝙱, 𝙲], 𝙱)] sage: duel_scores = [[33, 27], [25, 35], [42, 18]]  Individual scores: sage: scores = individual_duel_scores(duel_winners, duel_scores)  Extract the players from the scores: sage: players = sorted(scores) sage: players ['A', 'B', 'C']  Some scores: sage: scores[A][C] 25 sage: scores[A][A] 0  The scores as a matrix: sage: matrix(ZZ, 3, [scores[a][b] for a in players for b in players]) [ 0 33 25] [27 0 42] [35 18 0]  2021-01-19 14:52:16 -0600 commented question Tensor product of elements of non-free algebras The code a = SteenrodAlgebra(2).an_element() M = CombinatorialFreeModule(GF(2), ['s', 't', 'u']) s = M.basis()['s'] T = tensor([a, s])  works fine in Sage 8.8 but fails in Sage 8.9. This might have to do with the changes in Sage Trac ticket 25603. 2021-01-19 14:14:50 -0600 commented answer How to do low degree computation in a Free Algebra ? To format code blocks in questions, answers, or comments, separate them by a blank line above and below, and indent them by four spaces. Cay you edit your comment to do that? Or you could turn your comment into an answer, since it actually answers the question. 2021-01-19 14:11:58 -0600 answered a question Ranging the Z axis # Setting zmin and zmax to get a chosen z-range ## Missing functionality Inspecting the plot method of hyperplanes arrangement reveals this is not implemented in Sage yet. Searching Ask Sage and sage-support, we find that many 3d plotting functions similarly lack the ability to limit the z-range. Providing this missing functionality is now tracked at ## Workaround Some of the questions on this theme were answered with a suggestion to use implicit_plot3d. For hyperplanes, this works very well. We can therefore mimic the plot method of hyperplane arrangements and cook up a function which, if given only an x-range and a y-range, will simply call the plot method (which adapts the z-range to the full range of z values corresponding to the x-range and y-range), but if additionally given a z-range, will use implicit plots to limit the z-range. Here is such a function: def ha_plot(A, ranges, **opt): r""" Return a plot of this hyperplane arrangement. INPUT: - A -- a hyperplane arrangement - ranges -- a list consisting of an x-range, a y-range, and optionally a z-range """ from colorsys import hsv_to_rgb if len(ranges) == 2: return A.plot(ranges=ranges, **opt) elif len(ranges) != 3: raise ValueError("ranges should consist in two or three tuples") hh = [h.coefficients() for h in A.hyperplanes()] N = len(hh) HSV_tuples = [(i*1.0/N, 0.8, 0.9) for i in range(N)] cols = [hsv_to_rgb(*x) for x in HSV_tuples] ff = lambda a, b, c, d: lambda x, y, z: d + a*x + b*y + c*z ip = lambda f, col: implicit_plot3d(f, *ranges, color=col, **opt) hp = lambda a, b, c, d, col: ip(ff(a, b, c, d), col) return sum(hp(a, b, c, d, col) for (d, a, b, c), col in zip(hh, cols))  Usage is as follows. One can provide only an x-range and a y-range: sage: ha_plot(A, [(-1, 1), (-1, 1)]).show(aspect_ratio=(1, 1, 0.25))  and get the exact same picture as in the question. Or one can additionally provide a z-range: sage: ha_plot(A, [(-1, 1), (-1, 1), (-1, 1)]).show(aspect_ratio=1)  Optional extra arguments are passed to implicit_plot3d. This for instance allows to set the opacity: sage: ha_plot(A, [(-1, 1), (-1, 1), (-1, 1)], opacity=0.5).show(aspect_ratio=1)  The function could be improved to accept a choice of colours. Currently it uses the same default as the plot method of hyperplane arrangements. 2021-01-19 08:58:54 -0600 edited question Ranging the Z axis I am plotting arrangements of hyperplanes. Like this one: H3.<x, y, z> = HyperplaneArrangements(QQ) A = H3([(1, 2, 1), 0], [(-4, -3, -2), 0], [(3, 3, -2), 0], [(1, -4, 3), 0], [(-2, 2, 1), 0]) A.plot(ranges=[(-1, 1), (-1, 1)], aspect_ratio=(1, 1, 0.25))  The problem is I need the z axis to range from -1 to 1 rather than from -4 to 4. Is there a way to change that? 2021-01-19 08:55:04 -0600 commented question Is there a way to range the z axis? 2021-01-19 05:09:54 -0600 edited question Random polynomial of degree 1 Having defined • the finite field K of size 4 in a • a polynomial ring R1 over K one can ask for a random element using the random_element method. The random polynomial I obtained that way was of degree two: sage: K.<a> = FiniteField(4) sage: R1 = PolynomialRing(K, ['z%s' % p for p in range(1, 3)]) sage: R1.random_element() (a)*z1*z2 + (a + 1)*z2^2 + z1 + (a)*z2  If I only want a random linear polynomial from this ring, how can I get one? 2021-01-19 04:53:16 -0600 edited question Multiply polynomials from different rings Suppose I take a polynomial from$K[x]$, say$x^2 + 5x$, and another polynomial from$K[y]$, say$y^3$. I want to formally multiply them and get$x^2y^3 + 5xy^3$as the output. How can I do that? Note: K is the finite field of size 4 in a. My efforts: sage: K.<a> = FiniteField(4) sage: R_1 = PolynomialRing(K, ['z%s' % p for p in range(1, 3)]) sage: R_2 = PolynomialRing(K, ['x%s' % p for p in range(1, 3)]) sage: pp = R_1.random_element() sage: pp (a + 1)*z1^2 + (a)*z1*z2 + (a + 1)*z2 + (a) sage: qq = R_2.random_element() sage: qq x1*x2 + (a + 1)*x2^2 + x1 + 1  When I do pp*qq I get the following output unsupported operand parent(s) for *: 'Multivariate Polynomial Ring in z1, z2 over Finite Field in a of size 2^2' and 'Multivariate Polynomial Ring in x1, x2 over Finite Field in a of size 2^2'  2021-01-18 14:57:47 -0600 commented question "Affine diagonalization algorithm" in n-dimensions? The .mws file is really a text file, you can open it with a text editor and sort of make sense of it. 2021-01-18 14:54:25 -0600 edited question "Affine diagonalization algorithm" in n-dimensions? Does Sage have a general implementation of the "affine diagonalization algorithm" for n-dimensional vector spaces? I found some pseudo-code, see page 15 (in section 3 "affine diagonalization" which begins on page 12) of Searching online led me also to this answer from 2008 where someone does it with Maple: but unfortunately, I cannot open the .mws file (I get an error message concerning the version number). Thank you very much for the help. 2021-01-18 12:17:20 -0600 edited answer Centering the cells of a table and the use of tiny size Regarding the "\tiny doesn't work" part, you should add a $ before the LaTeX formula and you should either escape the \ with another \ or use a raw string:

sage: cond1 = ['$\\tiny{}{}{}{}$'.format(*w) for w in Ap]
sage: cond1 = [r'$\tiny{}{}{}{}$'.format(*w) for w in Ap]


Regarding the table, it needs a list of lists, not just a list, so use extra brackets; and you want headers for the columns, so they form a row of headers or "header row".

sage: t3 = table([total_score_de_Borda], header_row=cond3)
sage: t3
$A$   $B$   $C$   $D$
+-----+-----+-----+-----+
717   691   669   613

2021-01-18 12:08:18 -0600 edited question Centering the cells of a table and the use of tiny size

I change my question to clarify. In the following code (I have integrated the solution for \tiny):

First define a function:

def All_pref(cand=["A", "B", "C", "D"], code=1) :
ncand = len(cand)
Scand = sorted(Set(cand))
all_pref = Arrangements(Scand, ncand).list()
all_pref1 = [str(Word(x)) for x in all_pref]
if code == 1:
return ncand
if code == 2:
return all_pref1


Then use it:

cand = ["A", "B", "C", "D"]
ne1 = [18, 16, 14, 12, 11, 20, 19, 14, 16, 12, 2, 1, 0, 0, 20, 16, 13, 15, 11, 10, 9, 8, 7, 5]
ne = ne1
Ap = All_pref(cand, 2)
cond1 = ['$\\tiny${}{}{}{}$'.format(*w) for w in Ap] cond2 = [""] + [r'${}$'.format(*w) for w in cand] cond3 = [r'${}$'.format(*w) for w in cand] t1 = table([ne], header_row=cond1) rank = [[x.find(l) for x in Ap] for l in cand] score_de_Borda = [[abs(x.find(l) - len(cand)) for x in Ap] for l in cand] t2 = table(rank, header_row=cond1, header_column=cond2) total_score_de_Borda = [add([score_de_Borda[j][i]*ne[i] for i in range(24)]) for j in range(len(cand))] # t3 = table(total_score_de_Borda, header_column=cond3) show(cond3,total_score_de_Borda)  everything works fine apart from the header_column=cond3 (without it I have no problem) and as show(cond3, total_score_de_Borda) shows, the two elements have the same length. 2021-01-18 10:25:27 -0600 commented question I should store binary values in an array Why not open a new question: • click the "Edit" button under the question • copy the "follow-up" part • open a new question • paste the "follow-up" part there • remove it from here The new question can say "this is a follow up to Ask Sage question 55323" and you can add a comment here saying "See a follow-up question at Ask Sage question 55xxx". 2021-01-18 09:45:29 -0600 edited question How can I compute a fixed field over the p-adics I have problems to implement the following set up: I want to have a field$K = \mathbb{Q}_3$and an extension$L = \mathbb{Q}_3(\alpha)$over$K$where$f:=\min_K(\alpha) = x^4 - 3x^2 + 18$. This extension has degree$4$and ramification index$2$. Furthermore, let$F/K$be the unique unramified extension of$K$of degree$4$which is generated by a primitive$5$-th root of unity$\zeta_5$. Then one can show that$\varphi: \alpha \mapsto \frac{(2 \alpha^2 - 3)\sqrt{-\frac{2}{7}}}{\alpha},\zeta_5 \mapsto \zeta_5^3$is an element of the Galois group of$LF/K$. Now let$L' = (LF)^{\langle \varphi \rangle}$. This must be a quadratic and totally ramified extension of$K$. There are only two possibilities for that:$K(\sqrt{3})$or$K(\sqrt{-3})$. Question: How to determine whether$L' = K(\sqrt{3})$or$L' = K(\sqrt{-3})$(or equivalently,$\varphi(\sqrt{3}) = \sqrt{3}$or$\varphi(\sqrt{-3}) = \sqrt{-3}$)? Since I only have only superficial knowledge about Sage, I was not even able to set up the easy things like the field$L$properly. When I use K = Qp(3) R.<x> = ZZ[] f = x^4 - 3*x^2 + 18 L.<alpha> = K.extension(f)  I get an error that my polynomial$f$must be either unramified or Eisenstein (which of course does not exist since$L/K\$ is neither unramified nor totally ramified). Furthermore, I have no idea how to approach with my problem with Sage otherwise. And since computation by hand is pretty hard in this case (I already tried!), it would be nice to solve with problem here, so I can use it for similar computations in the future.