Dr. Patrick T

2083 Reputation

13 years, 259 days

maybe assuming and AllSolutions...

making assumptions about x may eliminate some candidate solutions, you may also  remove rootofs with the AllSolutions option in solve, maybe something like

```E[i, j] := Int(f[i](t)*u[j](t), t = -infinity .. infinity);
eqn[i, j] := u[j](x) = E[i, j];
k[i, j] := solve(eqn[i, j], x, AllSolutions) assuming x::positive;
```

subs...

also does the trick (but couldn't tell you which if any is to be commended)

```subs(Result,alpha[2]);
-1.536774530
```

manual method...

There must be a way to automate the following (with the use of the option discont=true), and someone else will certainly show how:

```restart;
with(plots):
with(plottools):

fdiscont((x^2-4*x+5)/(x-2),x=1..3);

[1.99972287449048180 .. 2.00074755004038485]

```

fdiscont helps you find the discontinuities of your function -- okay, so the asymptote is near 2 -- weird that Maple doesn't understand it's exactly 2.

Then I split the plot into two parts, so that Maple won't show a vertical line at the asymptote -- I then create my own asymptote by hand and give it the look I want (including dashes):

```p1:=plot((x^2-4*x+5)/(x-2), x=1..1.99):
p2:=plot((x^2-4*x+5)/(x-2), x=2.001..3):
p3:=curve([[2,-10], [2,10]], color=red, linestyle=dash, thickness=2):
display({p1,p2,p3}, view=[1..3,-10..10]);
```

maybe...

assuming you meant to type * instead of = in the formula for H:

```H:= (cos(n*L)*cosh(n*L))/(sin(n*L)-sinh(n*L)):
phi:=x->sin(n*x)-sinh(n*x)+H*(cos(n*x)-cosh(n*x)):
int(phi(x), x=0..L);
1/2 (4 exp(n L) sin(n L) - 4 exp(n L) sinh(n L)

- 2 cos(n L) sin(n L) exp(n L)

- cos(n L) cosh(n L) exp(2 n L)

+ 2 cos(n L) sinh(n L) exp(n L) + cos(n L) cosh(n L)

- exp(2 n L) sin(n L) + exp(2 n L) sinh(n L) - sin(n L)

+ sinh(n L) + 2 cos(n L) cosh(n L) sin(n L) exp(n L))

exp(-n L)/((sin(n L) - sinh(n L)) n)

```

exponential notation on the vertical axi...

loglogplot(10^x, x=1..2);

on my Windows/Maple13/Classic system I do get exponentials on the vertical axis.

You did not give us your code.

not very likely...

Explicit solutions are the exception rather than the rule. The following can be used to extract information from Maple's dsolve command.

```restart;
with(DEtools);
infolevel[dsolve] := 4;
```
```dsys1 := {diff(u(t), t) = -(3/2-sqrt(2))*u(t)+5*v(t), diff(v(t), t) = -5*u(t)-(3/2+sqrt(2))*v(t)-20*cos(10*t)*w(t), diff(w(t), t) = -.5-3*w(t)+20*cos(10*t)*v(t)};

dsolve(dsys1);

-> Solving each unknown as a function of the next ones using the order: [u(t), w(t), v(t)]
-> Calling odsolve with the ODE, diff(diff(diff(y(x),x),x),x) = -1/4*(-48000*cos(10*x)^2*y(x)*sin(10*x)+7980*cos(10*x)^2+3100*diff(y(x),x)*sin(10*x)-728*diff(diff(y(x),x),x)*cos(10*x)+600*diff(diff(y(x),x),x)*sin(10*x)-319200*cos(10*x)^3*y(x)-19291*y(x)*cos(10*x)+4800*cos(10*x)^3*diff(y(x),x)+9090*y(x)*sin(10*x)+32000*cos(10*x)^2*diff(y(x),x)*sin(10*x)-960000*cos(10*x)*y(x)*sin(10*x)^2-1989*diff(y(x),x)*cos(10*x)+16000*sin(10*x)^2+(1600*cos(10*x)*sin(10*x)-96000*cos(10*x)^2*y(x)*sin(10*x)+80*diff(diff(y(x),x),x)*sin(10*x)+240*diff(y(x),x)*sin(10*x)+2020*y(x)*sin(10*x)+3200*cos(10*x)^3*diff(y(x),x)+48*diff(diff(y(x),x),x)*cos(10*x)+274*diff(y(x),x)*cos(10*x)+606*y(x)*cos(10*x))*2^(1/2))/(3*cos(10*x)+20*sin(10*x)+2*2^(1/2)*cos(10*x)), y(x), singsol = none
Methods for third order ODEs:
--- Trying classification methods ---
trying high order exact linear fully integrable
trying differential order: 3; linear nonhomogeneous with symmetry [0,1]
trying high order linear exact nonhomogeneous
trying differential order: 3; missing the dependent variable
-> pFq: Equivalence to the 3F2 or one of its 3 confluent cases under 'a power @ Moebius'
trying a solution in terms of MeijerG functions
-> Try computing a Rational Normal Form for the given ODE...
<- unable to resolve the Equivalence to a Rational Normal Form
checking substitution methods
trying differential order: 3; missing the dependent variable
-> pFq: Equivalence to the 3F2 or one of its 3 confluent cases under 'a power @ Moebius'
trying a solution in terms of MeijerG functions
-> Try computing a Rational Normal Form for the given ODE...
<- unable to resolve the Equivalence to a Rational Normal Form
checking substitution methods
trying differential order: 3; missing the dependent variable
checking if the LODE is of Euler type
trying Louvillian solutions for 3rd order ODEs, imprimitive case
-> pFq: Equivalence to the 3F2 or one of its 3 confluent cases under 'a power @ Moebius'```

Maple then gives you an equivalent representation of your system in terms of one third-degree ODE in v(t). Run the above code and you can see that it's unwieldy.

I don't get the memory allocation error problem.

I am guessing that EU must be something like: Px*U(x)+Py*U(y), where U() is your utility function.

if your utility function is quadratic, since your constraint linear, you should/must be able to solve the problem in closed form.

if U(x)=1-x^2, say...

Axel's contourplot with numpoints=50^3...

Takes a short while to produce but looks pretty smooth:

```contourplot( max(min(-z*(1/(x^2+z^2)^(3/2)-1), 1000),-1000),
x=-0.3..0.3, z=-0.3..0.3, numpoints=50^3,axes=boxed);
```

explicit use of max and min...

The max min trick is pretty cool. Its usefulness is particularly apparent, in this example, for small values of numpoints:

```plot3d( max(min(-z*(1/(x^2+z^2)^(3/2)-1), 1000),-1000),x=-0.3..0.3,z=-0.3..0.3,
numpoints=50,axes=boxed,lightmodel=light2,orientation=[27, 54]);
```
```
```

versus

```plot3d( -z*(1/(x^2+z^2)^(3/2)-1), x=-0.3..0.3, z=-0.3..0.3,
view=[DEFAULT,DEFAULT,-1000..1000], numpoints=50,
axes=boxed,lightmodel=light2,orientation=[27, 54]);

```
```
```

Okay, I confess, I'm really practicing my image-posting skill.

oops, sorry!...

you're absolutely right, I copied the wrong line of code, apologies.

system explodes...

I think. You can see that by gradually raising the range of t from, say, 0..1 to 0..4, the following sums it up:

```plots[odeplot](sol,[t,Theta1(t)],0..2);
plots[odeplot](sol,[t,Theta2(t)],0..2);
```

solve...

in this case, solve does it:

y=(-4+4*d)*x^2+(4-4*d)*x-3*d+2;
solve(%, x);

numpoints=n, n>50...

one way to do it: numpoints

the default numpoints is 50, increase it for better effect:

contourplot(-z*(1/(x^2+z^2)^(3/2)-1),x=-3..3,z=-3..3,grid=[50,50],contours=10,numpoints=100);

is this what you're looking for?...

with(plots):
f := x -> x^(-2):
logplot(f(x), x=0..2);

 First 14 15 16 17 18 19 20 Last Page 16 of 24
﻿