Adri vanderMeer

## 19 Badges

18 years, 294 days
University of Twente (retired)
Enschede, Netherlands

My "website" consists of a Maple Manual in Dutch

## Use DEplot...

Convert the system to a system of first order equations:

`restart; with(DEtools):M:=1: L:=1: g:=9.81:sys := diff(x(t),t) = dx(t), diff(z(t),t) = dz(t),diff(dz(t),t)=-1/z(t)*diff(x(t),t)^2-1/z(t)*diff(z(t),t)^2-x(t)/z(t)*diff(dx(t),t),M*(-L^2/(x(t)*z(t)))*diff(dx(t),t)=M*g-M*(-1/z(t)*diff(x(t),t)^2-1/z(t)*diff(z(t),t)^2);`

Now use DEplot to plot z vs x:

```DEplot( {sys}, [x,z,dx,dz], t=0..2,
{[dx(0)=0,x(0)=0.99999999,dz(0)=0,z(0)=sqrt(L^2-0.99999999^2)]},
stepsize=0.001, scene=[x,z], x=-2..2, z=-.1..1.1, linecolor=blue );
```

DEplot.mw

## Determinants of minors...

Of course you can calculate the eigenvalues, and check if thesse are all positive. In this case these are rather complicated expressions. Easier - especially when your matrices become larger - is to calculate the determinant of submatrices. If these are all positive, the quadratic form is positive definite; if the signs of these determinants alternate, it is negative definite. So, for your example:

`restart;A:=Matrix([[8,4,2,1],[4,8,2,1],[2,2,8,1],[1,1,1,8]]):seq( LinearAlgebra:-Determinant(A[1..i,1..i]), i=1..4 );                        8, 48, 352, 2736`

All positive, so the quadratic form is positive definite.

## By making all the free parameters zero...

Your system is underdeterminated, zo there are many free parameters. A naive approach is to set all free parameters zero. I give a small example:

`restart; with(LinearAlgebra):A := RandomMatrix(4,9);b := <0,3,0,0>;c := LinearSolve(A,b);eval(c,indets(c)=~0);`

## output=Array...

Another possibility is:

`dsnx := dsolve(sysx, numeric, method = bvp[midrich], continuation = lambda,        maxmesh = 20000, initmesh = 20000, abserr = 10^(-2),        output = Array([seq(i,i=1..10,0.1)]) );res := dsnx[2,1]:res[[1,10,20,40,50,90],[1..3]];`

## Rather complicated......

I added de derivative of the third equation, and I make all variables functions of t:

`restart;a := miu1(t) = Diff(lambda1(t),t) + lambda3(t);b := miu2(t) = Diff(lambda2(t),t) + lambda1(t) - 4*lambda3(t);c := miu3(t) = Diff(lambda3(t),t) + lambda2(t) - 3*lambda3(t);c1 := Diff(lambda3(t), t, t) = -(Diff(lambda2(t), t)) + 3*(Diff(lambda3(t), t)) + Diff(miu3(t), t);d := miu(t) = 2*Diff(lambda3(t),t\$2) - 5*Diff(lambda3(t),t) + 7*lambda3(t);`

Now I use the first three equations to solve for the derivatives of lambda, and substitute rhe result in the last equation:

` s := solve( {a,b,c}, {Diff(lambda1(t),t),Diff(lambda2(t),t),Diff(lambda3(t),t)} ):A1 := subs(s,d);`

This contains the second derivative of lambda3, which is also in the third equation. But this equation, in turn, contains again derivatives of lambda. So I use the 2nd ant 3rd equation to get expressions for these derivatives in terms of lambda:

```A2 := solve( {b,c}, {Diff(lambda2(t),t),Diff(lambda3(t),t)} );
A3 := subs(A2,c1);
subs(A3,A1);
```

I don't understand your second question.

## Only two solutions...

Maple's answer is correct:

`restart;solve([-sin(x) = cos(x), 2*Pi > x and x > 0], x, AllSolutions = true);                     /      1            \                     { x = - - Pi + Pi _Z1 }                     \      4            / about(_Z1);Originally _Z1, renamed _Z1~:  is assumed to be: AndProp(integer,RealRange(1,2))`

This shows that _Z1=1 or _Z1=2.

## Syntax...

`DE := diff(y(x),x\$2) - 3*diff(y(x),x) - 4 = exp(x)*sin(4*x);dsolve(DE);`

Try it yourself the next time

## How about this?...

`normalcdf:= proc(x1,x2,mu,sigma)local t, dist:  dist:= t -> evalf(Statistics:-CDF(Normal(mu, sigma), t));  print(P(z<x1)=dist(x1));  print(P(z>x2)=1-dist(x2));  print(P(abs(z-mu)>sigma)=1-dist(x2)+dist(x1));  print(P(abs(z-mu)<sigma)=dist(x2)-dist(x1))end proc:normalcdf(2100,2200,2150,50);                   P(z < 2100) = 0.1586552540                   P(2200 < z) = 0.1586552540               P(50 < |z - 2150|) = 0.3173105080               P(|z - 2150| < 50) = 0.6826894920`

## use fsolve...

Notice that your function has two zeros at x=0 for all omega, so divide by x^2. To plot the other zero, depending on omega:

`restart;f := x + omega^2*x^3 - omega:pnts := [seq( [omega, fsolve(f)], omega=0..10,0.1 )]:plot( pnts );`

I suppose that you can write a Newton-Raphson procedure to use instead of fsolve if it is essential for your problem to use this method.

## check...

Your second question is to check the solution:

`restart:g:= x-> d*x^2+e*x+f:sol := solve({g(0)=4, g(3)=7, g(-2)=18}, {d,e,f}):G := subs( sol, eval(g) );     8  2   19      x -> - x  - -- x + 4     5      5       is( G(0)=4 and G(3)=7 and G(-2)=18 );                              true`

## Pi...

Replace 2sin(pix)  by 2*sin(Pi*x) : Pi with capital P, because pi is only a name for a variable.

## Straightforward...

If the list of x-values is called X, and the lists of y-values are called Y1,...Ym:

`L := [seq( seq([X[j],Y||j[i]], i=1..n), j=1..m )]:plots:-pointplot(L);`

## use floats...

replace lengte := 22 by lengte := 22.0

## Vector(field)?...

Your notation for f suggests that you mean a 3d vectorfield in cartesian coordinates. See ?VectorCalculus

Using this package, the unit vectors i, j, k are represented as ex, ey, ez.

## Parametrization...

You will get a nicer plot if you make a parametrization of the surface:

`plot3d( [13.2*sin(phi)*cos(theta), 24.3*sin(phi)*sin(theta), 22*cos(phi)],              phi=0..Pi, theta=0..2*Pi, scaling=constrained );`
 First 9 10 11 12 13 14 15 Last Page 11 of 27
﻿