Ask Your Question

calc314's profile - activity

2019-06-08 10:04:35 -0600 received badge  Popular Question (source)
2019-03-09 12:05:36 -0600 edited question TypeError: cannot solve differential equation

I am trying to use sagemath to solve a system of ordinary differential equation. However, I failed and got an TypeError (This returns the following error: TypeError: Computation failed since Maxima requested additional constraints; using the 'assume' command before evaluation may help). Here is my code for it. Could anyone help me with this? Thank you.

A=Matrix([[5/3, -4/3], [4/3, -5/3]])
B=Matrix([[2],[1]])
C=Matrix([-1,2])
D=Matrix([0])
t = var('t')
X=symb_mat(2,1,'X',t);
diffeqns = [eql == eqr \
            for (eql,eqr) in zip(*[diff(X,t).list(),(A*X+B*u).list()])]
desolve_system(diffeqns,X.list(),ics=[0,1,2],ivar=t)
2019-01-17 08:59:20 -0600 answered a question How can I solve the following (linear) differential equation?

sympy does give a solution, but it leaves you to solve the characteristic equation.

from sympy import *
y=Function('y')
dsolve(Eq(Derivative(y(x),x,x,x)-3*Derivative(y(x),x,x)+Derivative(y(x),x) -5*y(x), 0), y(x))
2019-01-09 12:51:05 -0600 received badge  Nice Answer (source)
2018-10-02 07:50:50 -0600 answered a question plotting autonomous differential equations

Here is another option for plotting. The streamline_plot does a nice job of plotting some sample curves along with arrows to give direction. I like to use the odeint solver since it works for systems so nicely. Here I use it by defining a system $ds/dt = 1, dx/dt = x^2-4x$. This gives the appropriate plot. Note that we have to tell Sage that the number 1 should be in the symbolic ring (SR) to get this to work.

var('x s t')
pp=streamline_plot(x^2-4*x,(t,0,5),(x,-4,5))
initialvalues=[-2,2]
for x0 in initialvalues:
    pp+=line(desolve_odeint([SR(1),x^2-4*x],[0,x0],srange(0,5,0.05),[s,x]),color='red',thickness=3)
pp.show(ymin=-4,ymax=5)
2018-08-31 20:30:08 -0600 received badge  Notable Question (source)
2018-07-18 05:54:33 -0600 received badge  Nice Answer (source)
2018-06-28 13:47:20 -0600 received badge  Nice Answer (source)
2018-06-11 12:57:11 -0600 received badge  Famous Question (source)
2018-05-08 20:13:06 -0600 answered a question One interact slider for multiple functions

I'm not sure exactly what you are aiming for, but in a sagews worksheet in CoCalc, you can do the following.

@interact
def go(f1=x+1,f2=x+2, a=slider(1,5,1,default = 3)):
    p1=plot(a*f1,(x,-3,3))
    p2=plot(f2+a,(x,-3,3))
    show(p1+p2)
2018-05-04 12:10:56 -0600 answered a question Solve a system

For the first solve question, you need to specify that you are solving for x and y.

solve([2*x+y==0,x-2*y==1],[x,y])

For the second question, you can do the following:

var('x y')
assume(x,'real')
assume(y,'real')
z=x+i*y
f=5*z-i*conjugate(z)^2+3+2*i
solve([f.real(),f.imag()],[x,y])
2018-04-17 10:20:56 -0600 received badge  Famous Question (source)
2018-04-11 11:31:56 -0600 received badge  Nice Question (source)
2018-04-09 23:34:00 -0600 received badge  Nice Answer (source)
2018-03-09 03:13:46 -0600 received badge  Famous Question (source)
2018-03-03 06:52:14 -0600 answered a question Is it possible to speed up loop iteration in Sage?

Using Cython can dramatically improve performance by converting the code to C and compiling it. For example,

%%cython
import random
def run(end): 
    uni = {}

    for count in range(0, end):
        i = random.randint(1, 303325737249669131)  
        if i in uni:
            uni[i] += 1
        else:
            uni[i] = 1
    return uni
2018-02-08 20:09:11 -0600 edited question How do get SVD inverse

I formed my matrix $A$ from, A = random_matrix(QQ, 6, algorithm='echelonizable', rank=6, upper_bound=9), I then changed the ring to RDF (A.change_ring(RDF)) , I further obtain my SVD matrix (U, S, V = B.SVD()). How do i find the inverse of A using SVD

2018-01-12 20:37:35 -0600 answered a question Finding extrema and zeros in lists

Another option is to use splines to approximate the curve. Here is some code inspired by https://stackoverflow.com/questions/3...

var('t y z')
P=desolve_system_rk4([z,-z^2*(-2+t/(1+t))-y],[y,z],ics=[0.1,10,2],ivar=t,end_points=100)
Q=[ [float(t1),float(z1)] for t1,y1,z1 in P]
x_points=([p[0] for p in Q])
y_points=([p[1] for p in Q])
tmp = interpolate.splrep(x_points, y_points)
f=lambda x: interpolate.splev(x, tmp)

You could then find a zero using:

find_root(f,12,20)

Also, you can compute the derivative using:

g=lambda x: interpolate.splev(x, tmp,der=1)

Then, find potential extrema using find_root.

2018-01-12 20:12:12 -0600 edited answer Finding extrema and zeros in lists

OK, I think I have a path to a solution. I have changed the differential equation so that it has zeros and this code prints the values on either side of where the function crosses zero. I can do something similar for the extrema. So I guess I am now asking if there is a function call that does this, or if there is a smarter way. Thanks, Wayne

var('t y z')

P=desolve_system_rk4([z,-z^2*(-2+t/(1+t))-y],[y,z],ics=[0.1,10,2],ivar=t,end_points=100)

Q=[ [t1,y1] for t1,y1,z1 in P]

line(Q).show()

print len(Q)

[a1,b1]=Q[0]

for i in range(1,len(Q)):

    [a2,b2]=Q[i]

    if(((b1<0.)and(b2>0.)) or ((b1>0.)and(b2<0.))):

        print a1,b1,a2,b2

    a1=a2

    b1=b2
2018-01-11 14:07:32 -0600 commented answer Plot/Show Png/jpg in sagemath worksheet

You can double-click the sell to unhide it.

2018-01-08 10:42:22 -0600 commented answer Plot/Show Png/jpg in sagemath worksheet

Try adding %hide before the %md.

2018-01-08 06:28:29 -0600 received badge  Nice Answer (source)
2018-01-07 12:16:35 -0600 answered a question Plot/Show Png/jpg in sagemath worksheet

I would create a markdown cell and link to the image using a markdown command. For example:

%md
![text description](example.jpg)
2017-12-30 06:49:28 -0600 commented answer solving nonlinear second order ordinary differential equations numerically

@eric_g's solution can handle this. Just replace $g(t)$ with $y(t)$. Here is an example with $f(t)=-2+\frac{t}{1+t}$.

var('t y z')
P=desolve_system_rk4([z,-z^2*(-2+t/(1+t))-y],[y,z],ics=[0.1,10,2],ivar=t,end_points=100)
Q=[ [t1,y1] for t1,y1,z1 in P]
line(Q).show()
2017-12-27 07:31:22 -0600 edited question Why does the code not work?

Hi, I wanted to use a program from an older question from myself again and found the following strange thing:

The following code works fine:

m=3

le_relations = lambda P: [(a,b) if P.le(a,b) else (b,a) for a,b in P.comparability_graph().edges(labels=None)]

cover_relations = lambda P: [(a,b) if P.le(a,b) else (b,a) for a,b in P.hasse_diagram().edges(labels=None)]

format_pt = lambda k: "'x{0}'".format(k+1)

format_relations = lambda relations: [[format_pt(a),format_pt(b)] for a,b in relations]

format_pts = lambda P: [format_pt(a) for a in P]

format_poset = lambda P: [format_pts(P), format_relations(cover_relations(P)), format_relations(le_relations(P))]

what_you_want = lambda n: [format_poset(P) for P in [posets.BooleanLattice(m)]]

print(what_you_want(m))

The following code does not work:

m=3

le_relations = lambda P: [(a,b) if P.le(a,b) else (b,a) for a,b in P.comparability_graph().edges(labels=None)]

cover_relations = lambda P: [(a,b) if P.le(a,b) else (b,a) for a,b in P.hasse_diagram().edges(labels=None)]

format_pt = lambda k: "'x{0}'".format(k+1)

format_relations = lambda relations: [[format_pt(a),format_pt(b)] for a,b in relations]

format_pts = lambda P: [format_pt(a) for a in P]

format_poset = lambda P: [format_pts(P), format_relations(cover_relations(P)), format_relations(le_relations(P))]

what_you_want = lambda n: [format_poset(P) for P in posets.YoungsLattice(m)]

print(what_you_want(m))

The only difference between two codes is that in the second code I have

what_you_want = lambda n: [format_poset(P) for P in posets.YoungsLattice(m)]

instead of

what_you_want = lambda n: [format_poset(P) for P in [posets.BooleanLattice(m)]].

Entering the second code gives an error. I do not really understand why the first works and the second does not, so any help is welcome how to fix the error.

2017-12-17 14:17:58 -0600 answered a question how to run maxima code in Sage?

Try placing the magic %maxima at the beginning of a Sage cell. Then, you can execute Maxima commands.

2017-12-17 14:15:11 -0600 edited question how to run maxima code in Sage?

I have the following code in maxima to calculate the laplacian of a function in parabolic coordinates.

assume(r≥0)$

assume(theta≥0,theta≤2*π)$

load(vect)$

derivabbrev:true$

scalefactors(parabolic)$

declare(f,scalar)$

depends(f,rest(parabolic))$

ev(express(laplacian(f)),diff,expand,factor);

What is the proper way to run those Maxima commands in SageMath?

Thanks,

Daniel

2017-12-01 12:33:52 -0600 answered a question What is useful for to upgrade a project or an account?
  1. The trial servers tend to be heavily used and so will run a bit slow. The member servers are significantly faster!
  2. Yes, the network access is internet access. The membership allows network access -- you can download files directly into your account or use ssh to move data between your local computer and CoCalc.
2017-10-21 22:12:28 -0600 answered a question how to print the printout not too long on Sagecell

How about:

num_rows=4
for i in range(num_rows):
    print [myList[i*64/num_rows+j] for j in range(64/num_rows)]
2017-10-15 13:22:56 -0600 edited question plot operation error

I'm trying to plot the following.

def myfn2(x): if x<0: return 1 else: return -1 plot(myfn2(x),x,-3,3,figsize=3,color="red")

The graph is only displayed as -1. Why?

2017-10-08 14:20:14 -0600 received badge  Notable Question (source)
2017-10-07 13:51:29 -0600 commented question Get the constant value of an equation

You can define f(x,y,z)=2*x + 2*y + 2*z - 51/4 and then find f(0,0,0).

2017-10-07 13:45:17 -0600 answered a question Solving a system of equations -- small known number hinders the solution

Here is an approach using Sympy.

from sympy.solvers import solve
from sympy import Symbol
x=Symbol('x')
y=Symbol('y')
z=Symbol('z')

K=4.48e-13
xi=0.75
yi=0
zi=0

eq1=y**2*z/x**2-K
eq2=xi+yi-(x+y)
eq3=2*xi+yi+2*zi-(2*x+y+2*z)

print eq1,eq2,eq3

ans=solve([eq1,eq2,eq3],x,y,z)
print ans
2017-10-07 13:12:48 -0600 edited question Solving a system of equations -- small known number hinders the solution

Hello, I'm trying to solve a system of equations in the treatment of an equilibrium system. I don't get any answers when I use the code below as is.

If I define K = 4.48e-13:

var('x,y,z')
K=4.48e-13
xi=0.75
yi=0
zi=0
eq1=K==y^2*z/x^2
eq2=xi+yi==x+y
eq3=2*xi+yi+2*zi==2*x+y+2*z
solns = solve([eq1,eq2,eq3],x,y,z,solution_dict=True)
[[s[x].n(30), s[y].n(30), s[z].n(30)] for s in solns]

I get:

[]

I do get numbers when I change the known value of "K" in the 2nd line to 4.48e-10 or bigger.

var('x,y,z')
K=4.48e-10
xi=0.75
yi=0
zi=0
eq1=K==y^2*z/x^2
eq2=xi+yi==x+y
eq3=2*xi+yi+2*zi==2*x+y+2*z
solns = solve([eq1,eq2,eq3],x,y,z,solution_dict=True)
[[s[x].n(30), s[y].n(30), s[z].n(30)] for s in solns]

[[0.75039762 - 0.00068968040I, -0.00039762382 + 0.00068968023I, -0.00019881201 + 0.00034484028I], [0.75039762 + 0.00068968040I, -0.00039762382 - 0.00068968023I, -0.00019881201 - 0.00034484028I], [0.74920475, 0.00079524857, 0.00039762445]]

I did not test past 1e-10. I tried increasing the number of bits to 100 for all 3 unknowns with no luck. Is there any way to solve the problem using SageMath using the smaller number? Thanks!

2017-09-28 12:16:25 -0600 commented answer Displaying images with matplotlib

This works in a .sagews file in CoCalc. It does not work interactively in a Jupiter notebook in CoCalc, but I understand that interacts don't work in that notebook in CoCalc right now.

2017-09-23 17:39:59 -0600 edited question Output too long in SAGE

I have a program whos output is too long to be displayed in the terminal where I run sage. Is there a possiblity that the output of SAGE is automatically "copied" (as if I mark it and use ctrl+c) and then I can paste it into a text file? Or is there another way that the result can be automatically be saved into a text file or a way to make the terminal allow larger text before the text gets cut off?

Here the code Im talking about:

le_relations = lambda P: [(a,b) if P.le(a,b) else (b,a) for a,b in P.comparability_graph().edges(labels=None)]
cover_relations = lambda P: [(a,b) if P.le(a,b) else (b,a) for a,b in P.hasse_diagram().edges(labels=None)]

format_pt = lambda k: "'x{0}'".format(k+1)
format_relations = lambda relations: [[format_pt(a),format_pt(b)] for a,b in relations]
format_pts = lambda P: [format_pt(a) for a in P]
format_poset = lambda P: [format_pts(P), format_relations(cover_relations(P)), format_relations(le_relations(P))]

what_you_want = lambda n: [format_poset(P) for P in posets(n) if P.is_lattice()]

print(what_you_want(8))

It takes about 20 minits for my computer to finish the calculation while the calculation seems to never finish in the online sage cell https://sagecell.sagemath.org/.

The code is from the thread https://ask.sagemath.org/question/388... )

edit: Another option would be to divide the set [format_poset(P) for P in posets(n) if P.is_lattice] into smaller pieces (lets say 5 pieces) and apply the code to those smaller pieces. The set contains 222 elements. Is there a way to apply the cod to the first 40 of those elements first and then I do it seperately for the next 40 and so on till I have all. Maybe it is possible to filter instead of posets(n) (for n=8) over just the first 40 posets of posets(n) and then the next 40 and so on.

2017-09-16 16:38:59 -0600 edited question solution to logistic differential y'=r*y*(1-y)/b

How to solve this logistic differential in sage?

y=function('y')(x)
(b,r,y0)=(30,0.8,1.5)
de=diff(y,x)==r*y*((1-y)/b) # logistic equation
soln=desolve(de,y,ics=[0,y0])

mathematica solution is: y(x)=(b*exp(r*x)*y0)/(exp(r*x)*y0-y0+b)

Other questions/answers found on the logistic topic were not helpful for this problem

2017-09-14 07:43:53 -0600 marked best answer How to color a 3d plot by z-level?

I'd like to color a 3d plot based on z-level. I think this is easy to do in Maple or Mathematica, but I've been searching the web for help on doing this in Sage and can't find anything to help with plot3d or implicit_plot3d.

Here's the implicit_plot3d I'm using.

var('x,y,z')    
implicit_plot3d(x^2-y^2*z == 0,(x,-4,4),(y,-4,4),(z,-4,4)).show(mesh=True)

Also, are there color maps in Sage that produce plots with colors and lighting similar to the default in Mathematica?

2017-09-12 10:59:02 -0600 answered a question plot sage Graphics() from jupyter notebok in CoCalc

Here is something that works, although I don't view it as perfect.

from sage.all import *
from IPython.display import Image
var('x')
p=plot(sin(x),(x,-2*pi,2*pi))
p.save('tmp.png')
Image("tmp.png")
2017-09-08 08:52:24 -0600 answered a question hi , what wrong with this code

Sage wants to put two tickmarks on the graph and is asking that you extend the x-range so that it can do so. Try changing the code to:

show(plot(sin(a*x+b), (x,0,a*pi),ticks=pi/2*a,tick_formatter='latex',gridlines=True))
2017-08-02 16:50:45 -0600 received badge  Popular Question (source)
2017-07-27 09:22:42 -0600 received badge  Nice Answer (source)
2017-07-13 20:41:10 -0600 received badge  Good Question (source)
2017-07-10 12:25:14 -0600 answered a question Collecting coefficients of derivatives

Try this:

p=a1*a2*a3*x*f(x) + a1*a2*x^2*diff(f(x), x) + a1*a3*x^2*diff(f(x), x) + a2*a3*x^2*diff(f(x), x)
p.collect(diff(f(x),x))