calc314's profile - activity

Is there a way to set the kind of axes shown in a 3d plot? Seems like you should be able to specify boxed, axes, etc. in the show command.

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)

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))

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)

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)

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])

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?

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

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

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.

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

You can double-click the sell to unhide it.

Try adding %hide before the %md.

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

%md
![text description](example.jpg)

@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()

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.

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

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

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.

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)]

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?

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

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

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!